Changeset 876


Ignore:
Timestamp:
05/28/19 22:44:58 (5 years ago)
Author:
jisesh
Message:

devel: moved DYSL into compute_caldyn_slow_NH.F90 and compute_caldyn_Coriolis.F90

Location:
codes/icosagcm/devel
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • codes/icosagcm/devel/Python/src/kernels_caldyn_NH.jin

    r821 r876  
    250250END_BLOCK 
    251251 
    252 KERNEL(caldyn_slow_NH) 
    253   FORALL_CELLS_EXT('1', 'llm+1') 
    254      ON_PRIMAL 
    255         CST_IF(IS_INNER_INTERFACE, w_il(CELL) = 2.*W(CELL)/(rhodz(KDOWN(CELL))+rhodz(KUP(CELL))) ) 
    256         CST_IF(IS_BOTTOM_LEVEL,    w_il(CELL) = 2.*W(CELL)/rhodz(KUP(CELL)) ) 
    257         CST_IF(IS_TOP_INTERFACE,   w_il(CELL) = 2.*W(CELL)/rhodz(KDOWN(CELL)) ) 
    258      END_BLOCK 
    259   END_BLOCK 
    260   FORALL_CELLS_EXT('1', 'llm+1') 
    261      ON_EDGES 
    262        ! compute DePhi, v_el, G_el, F_el 
    263        ! v_el, W2_el and therefore G_el incorporate metric factor le_de 
    264        ! while DePhil, W_el and F_el do not 
    265        W_el         = .5*( W(CELL2)+W(CELL1) ) 
    266        DePhil(EDGE) = SIGN*(Phi(CELL2)-Phi(CELL1)) 
    267        F_el(EDGE)   = DePhil(EDGE)*W_el 
    268        W2_el        = .5*LE_DE * ( W(CELL1)*w_il(CELL1) + W(CELL2)*w_il(CELL2) ) 
    269        v_el(EDGE)   = .5*LE_DE*(u(KUP(EDGE))+u(KDOWN(EDGE))) ! checked 
    270        G_el(EDGE)   = v_el(EDGE)*W_el - DePhil(EDGE)*W2_el 
    271      END_BLOCK 
    272   END_BLOCK 
    273  
    274   FORALL_CELLS_EXT('1', 'llm+1') 
    275     ! compute GradPhi2, dPhi, dW 
    276     ON_PRIMAL 
    277        gPhi2=0. 
    278        dP=0. 
    279        divG=0 
    280        FORALL_EDGES        
    281          gPhi2 = gPhi2 + LE_DE*DePhil(EDGE)**2 
    282          dP = dP + LE_DE*DePhil(EDGE)*v_el(EDGE) 
    283          divG = divG + SIGN*G_el(EDGE) ! -div(G_el), G_el already has le_de 
    284        END_BLOCK 
    285        gradPhi2(CELL) = 1./(2.*AI) * gPhi2 
    286        dPhi(CELL) = gradPhi2(CELL)*w_il(CELL) - 1./(2.*AI)*dP 
    287        dW(CELL) = (-1./AI)*divG 
    288     END_BLOCK 
    289   END_BLOCK 
    290  
    291   ! We need a barrier here because we compute gradPhi2, F_el and w_il above and do a vertical average below 
    292   BARRIER 
    293  
    294   FORALL_CELLS_EXT() 
    295     ! Compute berni at scalar points 
    296     ON_PRIMAL 
    297        u2=0. 
    298        FORALL_EDGES        
    299          u2 = u2 + LE_DE*u(EDGE)**2 
    300        END_BLOCK 
    301        berni(CELL) = 1./(4.*AI) * u2 - .25*( gradPhi2(CELL)*w_il(CELL)**2 + gradPhi2(UP(CELL))*w_il(UP(CELL))**2 ) 
    302     END_BLOCK 
    303   END_BLOCK 
    304  
    305   FORALL_CELLS_EXT() 
    306     ON_EDGES 
    307     ! Compute mass flux and grad(berni) 
    308       uu = .5*(rhodz(CELL1)+rhodz(CELL2))*u(EDGE) - .5*( F_el(EDGE)+F_el(UP(EDGE)) ) 
    309       hflux(EDGE) = LE_DE*uu 
    310       du(EDGE) = SIGN*(berni(CELL1)-berni(CELL2)) 
    311     END_BLOCK 
    312   END_BLOCK 
    313  
    314 END_BLOCK 
  • codes/icosagcm/devel/Python/src/kernels_caldyn_hevi.jin

    r858 r876  
    127127END_BLOCK 
    128128 
    129 KERNEL(coriolis) 
    130 ! 
    131   DO iq=1,nqdyn 
    132     FORALL_CELLS_EXT() 
    133       ON_EDGES 
    134         Ftheta(EDGE) = .5*(theta(CELL1,iq)+theta(CELL2,iq))*hflux(EDGE) 
    135       END_BLOCK 
    136     END_BLOCK 
    137     FORALL_CELLS() 
    138       ON_PRIMAL 
    139         divF=0. 
    140         FORALL_EDGES 
    141           divF = divF + Ftheta(EDGE)*SIGN 
    142         END_BLOCK 
    143         dtheta_rhodz(CELL,iq) = -divF / AI 
    144       END_BLOCK 
    145     END_BLOCK 
    146   END DO ! iq 
    147 ! 
    148   FORALL_CELLS() 
    149     ON_PRIMAL 
    150       divF=0. 
    151       FORALL_EDGES 
    152         divF = divF + hflux(EDGE)*SIGN 
    153       END_BLOCK 
    154       convm(CELL) = -divF / AI 
    155     END_BLOCK 
    156   END_BLOCK 
    157 ! 
    158   FORALL_CELLS()  
    159     ON_EDGES 
    160       du_trisk=0. 
    161       FORALL_TRISK 
    162         du_trisk = du_trisk + WEE*hflux(EDGE_TRISK)*(qu(EDGE)+qu(EDGE_TRISK)) 
    163       END_BLOCK 
    164       du(EDGE) = du(EDGE) + .5*du_trisk 
    165     END_BLOCK 
    166   END_BLOCK 
    167  
    168 END_BLOCK 
  • codes/icosagcm/devel/src/dynamics/compute_caldyn_Coriolis.F90

    r844 r876  
    44  PRIVATE 
    55 
     6#include "../unstructured/unstructured.h90" 
     7 
    68  PUBLIC :: compute_caldyn_Coriolis 
    79 
    810CONTAINS 
     11 
     12#ifdef BEGIN_DYSL 
     13 
     14KERNEL(coriolis) 
     15! 
     16  DO iq=1,nqdyn 
     17    FORALL_CELLS_EXT() 
     18      ON_EDGES 
     19        Ftheta(EDGE) = .5*(theta(CELL1,iq)+theta(CELL2,iq))*hflux(EDGE) 
     20      END_BLOCK 
     21    END_BLOCK 
     22    FORALL_CELLS() 
     23      ON_PRIMAL 
     24        divF=0. 
     25        FORALL_EDGES 
     26          divF = divF + Ftheta(EDGE)*SIGN 
     27        END_BLOCK 
     28        dtheta_rhodz(CELL,iq) = -divF / AI 
     29      END_BLOCK 
     30    END_BLOCK 
     31  END DO ! iq 
     32! 
     33  FORALL_CELLS() 
     34    ON_PRIMAL 
     35      divF=0. 
     36      FORALL_EDGES 
     37        divF = divF + hflux(EDGE)*SIGN 
     38      END_BLOCK 
     39      convm(CELL) = -divF / AI 
     40    END_BLOCK 
     41  END_BLOCK 
     42! 
     43  FORALL_CELLS() 
     44    ON_EDGES 
     45      du_trisk=0. 
     46      FORALL_TRISK 
     47        du_trisk = du_trisk + WEE*hflux(EDGE_TRISK)*(qu(EDGE)+qu(EDGE_TRISK)) 
     48      END_BLOCK 
     49      du(EDGE) = du(EDGE) + .5*du_trisk 
     50    END_BLOCK 
     51  END_BLOCK 
     52 
     53END_BLOCK 
     54 
     55#endif END_DYSL 
     56 
     57  SUBROUTINE compute_coriolis_unst(hflux,theta,qu,Ftheta, convm,dtheta_rhodz,du) 
     58    USE ISO_C_BINDING, only : C_DOUBLE, C_FLOAT 
     59    USE grid_param, ONLY : nqdyn 
     60    USE data_unstructured_mod, ONLY : id_coriolis, primal_num, dual_num, edge_num, & 
     61          left, right,primal_deg,primal_edge,primal_ne,trisk_deg,wee,trisk,Ai, & 
     62          enter_trace, exit_trace 
     63    FIELD_U     :: hflux, Ftheta, qu, du 
     64    FIELD_MASS  :: convm 
     65    FIELD_THETA :: theta, dtheta_rhodz 
     66    DECLARE_INDICES 
     67    DECLARE_EDGES 
     68    NUM :: divF, du_trisk 
     69    START_TRACE(id_coriolis, 3,4,0) ! primal, dual, edge 
     70#include "../kernels_unst/coriolis.k90" 
     71    STOP_TRACE 
     72  END SUBROUTINE compute_coriolis_unst 
    973 
    1074  SUBROUTINE compute_caldyn_Coriolis(hflux,theta,qu, convm,dtheta_rhodz,du) 
  • codes/icosagcm/devel/src/dynamics/compute_caldyn_slow_NH.F90

    r851 r876  
    44  PRIVATE 
    55 
     6#include "../unstructured/unstructured.h90" 
     7 
    68  PUBLIC :: compute_caldyn_slow_NH 
    79 
    810CONTAINS 
     11 
     12#ifdef BEGIN_DYSL 
     13 
     14KERNEL(caldyn_slow_NH) 
     15  FORALL_CELLS_EXT('1', 'llm+1') 
     16     ON_PRIMAL 
     17        CST_IF(IS_INNER_INTERFACE, w_il(CELL) = 2.*W(CELL)/(rhodz(KDOWN(CELL))+rhodz(KUP(CELL))) ) 
     18        CST_IF(IS_BOTTOM_LEVEL,    w_il(CELL) = 2.*W(CELL)/rhodz(KUP(CELL)) ) 
     19        CST_IF(IS_TOP_INTERFACE,   w_il(CELL) = 2.*W(CELL)/rhodz(KDOWN(CELL)) ) 
     20     END_BLOCK 
     21  END_BLOCK 
     22  FORALL_CELLS_EXT('1', 'llm+1') 
     23     ON_EDGES 
     24       ! compute DePhi, v_el, G_el, F_el 
     25       ! v_el, W2_el and therefore G_el incorporate metric factor le_de 
     26       ! while DePhil, W_el and F_el do not  
     27       W_el         = .5*( W(CELL2)+W(CELL1) ) 
     28       DePhil(EDGE) = SIGN*(Phi(CELL2)-Phi(CELL1)) 
     29       F_el(EDGE)   = DePhil(EDGE)*W_el 
     30       W2_el        = .5*LE_DE * ( W(CELL1)*w_il(CELL1) + W(CELL2)*w_il(CELL2) ) 
     31       v_el(EDGE)   = .5*LE_DE*(u(KUP(EDGE))+u(KDOWN(EDGE))) ! checked 
     32       G_el(EDGE)   = v_el(EDGE)*W_el - DePhil(EDGE)*W2_el 
     33     END_BLOCK 
     34  END_BLOCK 
     35 
     36  FORALL_CELLS_EXT('1', 'llm+1') 
     37    ! compute GradPhi2, dPhi, dW 
     38    ON_PRIMAL 
     39       gPhi2=0. 
     40       dP=0. 
     41       divG=0 
     42       FORALL_EDGES     
     43         gPhi2 = gPhi2 + LE_DE*DePhil(EDGE)**2 
     44         dP = dP + LE_DE*DePhil(EDGE)*v_el(EDGE) 
     45         divG = divG + SIGN*G_el(EDGE) ! -div(G_el), G_el already has le_de 
     46       END_BLOCK 
     47       gradPhi2(CELL) = 1./(2.*AI) * gPhi2 
     48       dPhi(CELL) = gradPhi2(CELL)*w_il(CELL) - 1./(2.*AI)*dP 
     49       dW(CELL) = (-1./AI)*divG 
     50    END_BLOCK 
     51  END_BLOCK 
     52 
     53  ! We need a barrier here because we compute gradPhi2, F_el and w_il above and do a vertical average below 
     54  BARRIER 
     55 
     56  FORALL_CELLS_EXT() 
     57    ! Compute berni at scalar points 
     58    ON_PRIMAL 
     59       u2=0. 
     60       FORALL_EDGES     
     61         u2 = u2 + LE_DE*u(EDGE)**2 
     62       END_BLOCK 
     63       berni(CELL) = 1./(4.*AI) * u2 - .25*( gradPhi2(CELL)*w_il(CELL)**2 + gradPhi2(UP(CELL))*w_il(UP(CELL))**2 ) 
     64    END_BLOCK 
     65  END_BLOCK 
     66 
     67  FORALL_CELLS_EXT() 
     68    ON_EDGES 
     69    ! Compute mass flux and grad(berni) 
     70      uu = .5*(rhodz(CELL1)+rhodz(CELL2))*u(EDGE) - .5*( F_el(EDGE)+F_el(UP(EDGE)) ) 
     71      hflux(EDGE) = LE_DE*uu 
     72      du(EDGE) = SIGN*(berni(CELL1)-berni(CELL2)) 
     73    END_BLOCK 
     74  END_BLOCK 
     75 
     76END_BLOCK 
     77 
     78#endif END_DYSL 
     79 
     80  SUBROUTINE compute_caldyn_vert_NH_unst(mass,geopot,W,wflux, eta_dot,wcov,W_etadot, du,dPhi,dW) 
     81    USE ISO_C_BINDING, only : C_DOUBLE, C_FLOAT 
     82    USE data_unstructured_mod, ONLY : left,right,edge_num,primal_num,dual_num,id_vert_NH, & 
     83        enter_trace, exit_trace 
     84    FIELD_MASS   :: mass, eta_dot, wcov, W_etadot ! IN, BUF*3 
     85    FIELD_GEOPOT :: geopot,W,wflux,dPhi,dW ! IN*3, INOUT*2 
     86    FIELD_U      :: du ! INOUT 
     87    DECLARE_INDICES 
     88    NUM :: w_ij, wflux_ij 
     89    START_TRACE(id_vert_NH, 6,0,1) 
     90#include "../kernels_unst/caldyn_vert_NH.k90" 
     91    STOP_TRACE 
     92  END SUBROUTINE compute_caldyn_vert_NH_unst 
    993 
    1094  SUBROUTINE compute_caldyn_slow_NH(u,rhodz,Phi,W, F_el,gradPhi2,w_il, hflux,du,dPhi,dW) 
Note: See TracChangeset for help on using the changeset viewer.