Changeset 913 for codes/icosagcm/devel/src/diagnostics/compute_pression.F90

Timestamp:

06/17/19 17:38:58 (5 years ago)

Author:

dubos

Message:

devel : compute_pression for unstructured mesh

File:

• codes/icosagcm/devel/src/diagnostics/compute_pression.F90 (moved) (moved from codes/icosagcm/devel/src/diagnostics/pression.f90) (3 diffs)

codes/icosagcm/devel/src/diagnostics/compute_pression.F90

@@ -1 @@
-MODULE pression_mod
+MODULE compute_pression_mod
+  USE compute_diagnostics_mod
+  USE icosa
+  USE omp_para
+  USE disvert_mod, ONLY : ap, bp, ap_bp_present
+  IMPLICIT NONE
+  PRIVATE
+
+#include "../unstructured/unstructured.h90"
+
+  PUBLIC :: pression, compute_pression_hex, compute_pression_unst, &
+       pression_mid, compute_pression_mid_hex, compute_pression_mid_unst
 
 CONTAINS
 
+#ifdef BEGIN_DYSL
+
+{%- macro compute_pression(llmax)%} 
+{%- set inner_loop=caller() %}
+{%- set llmax="'%s'"%llmax %}
+IF(ap_bp_present) THEN
+   IF(offset>0) THEN
+     FORALL_CELLS_EXT('1',{{ llmax }})
+       ON_PRIMAL
+         {{ inner_loop }}
+       END_BLOCK
+     END_BLOCK
+  ELSE
+     FORALL_CELLS('1',{{ llmax }})
+       ON_PRIMAL
+         {{ inner_loop }}
+       END_BLOCK
+     END_BLOCK
+  END IF
+END IF
+{%- endmacro %}
+
+KERNEL(compute_pression)
+{% call compute_pression('llm+1') %}
+p(CELL) = AP(CELL) + BP(CELL) * ps(HIDX(CELL))  
+{% endcall %}
+END_BLOCK
+
+KERNEL(compute_pmid)
+{% call compute_pression('llm') %}
+pmid(CELL) = .5*(AP(CELL)+AP(UP(CELL)) + (BP(CELL)+BP(UP(CELL))) * ps(HIDX(CELL)) )  
+{% endcall %}
+END_BLOCK
+
+#endif END_DYSL
+
   SUBROUTINE pression(f_ps,f_p)
-  USE icosa
-  IMPLICIT NONE
     TYPE(t_field), POINTER :: f_ps(:)
     TYPE(t_field), POINTER :: f_p(:)
@@ -27 +72 @@
 
   SUBROUTINE pression_mid(f_ps,f_pmid)
-  USE icosa
-  IMPLICIT NONE
     TYPE(t_field), POINTER :: f_ps(:)
     TYPE(t_field), POINTER :: f_pmid(:)
@@ -49 +92 @@
   END SUBROUTINE pression_mid
 
-  SUBROUTINE compute_pression(ps,p,offset)
-  USE icosa
-  USE disvert_mod
-  USE omp_para
-  IMPLICIT NONE
+!------------- hexagonal-mesh compute kernels --------
+
+#define AP(ij,l) ap(l)
+#define BP(ij,l) bp(l)
+
+  SUBROUTINE compute_pression_hex(ps,p,offset)
     REAL(rstd),INTENT(IN) :: ps(iim*jjm)
     REAL(rstd),INTENT(OUT) :: p(iim*jjm,llm+1)
     INTEGER,INTENT(IN) :: offset
-    INTEGER :: i,j,ij,l
-
-    IF(ap_bp_present) THEN
-      DO    l    = ll_begin, ll_endp1
-!      DO    l    = 1, llm + 1
-        DO j=jj_begin-offset,jj_end+offset
-          DO i=ii_begin-offset,ii_end+offset
-            ij=(j-1)*iim+i
-            p(ij,l) = ap(l) + bp(l) * ps(ij)
-          ENDDO
-        ENDDO
-      ENDDO
-    END IF
-
-  END SUBROUTINE compute_pression
+    INTEGER :: ij,l
+#include "../kernels_hex/compute_pression.k90"
+  END SUBROUTINE compute_pression_hex
   
-  SUBROUTINE compute_pression_mid(ps,pmid,offset)
-  USE icosa
-  USE disvert_mod
-  USE omp_para
-  IMPLICIT NONE
+  SUBROUTINE compute_pression_mid_hex(ps,pmid,offset)
     REAL(rstd),INTENT(IN) :: ps(iim*jjm)
     REAL(rstd),INTENT(OUT) :: pmid(iim*jjm,llm)
     INTEGER,INTENT(IN) :: offset
-    INTEGER :: i,j,ij,l
+    INTEGER :: ij,l
+#include "../kernels_hex/compute_pmid.k90"
+  END SUBROUTINE compute_pression_mid_hex
 
-    IF(ap_bp_present) THEN
-      DO    l    = ll_begin, ll_end
-        DO j=jj_begin-offset,jj_end+offset
-          DO i=ii_begin-offset,ii_end+offset
-            ij=(j-1)*iim+i
-            pmid(ij,l) = 0.5*(ap(l)+ap(l+1) + (bp(l)+bp(l+1)) * ps(ij))
-          ENDDO
-        ENDDO
-      ENDDO
-    END IF
+#undef AP
+#undef BP
 
-  END SUBROUTINE compute_pression_mid
+!----------- unstructured-mesh compute kernels --------
 
-END MODULE pression_mod
+#define AP(l,ij) ap(l)
+#define BP(l,ij) bp(l)
+  
+  SUBROUTINE compute_pression_unst(ps, p, offset)
+    FIELD_PS,     INTENT(IN)  :: ps
+    FIELD_GEOPOT, INTENT(OUT) :: p
+    INTEGER,      INTENT(IN)  :: offset
+    DECLARE_INDICES
+#include "../kernels_unst/compute_pression.k90"
+  END SUBROUTINE compute_pression_unst
+
+  SUBROUTINE compute_pression_mid_unst(ps, pmid, offset)
+    FIELD_PS,   INTENT(IN)  :: ps
+    FIELD_MASS, INTENT(OUT) :: pmid
+    INTEGER,    INTENT(IN)  :: offset
+    DECLARE_INDICES
+#include "../kernels_unst/compute_pmid.k90"
+  END SUBROUTINE compute_pression_mid_unst
+
+#undef AP
+#undef BP
+
+END MODULE compute_pression_mod