1 2! KGEN-generated Fortran source file 3! 4! Filename : prim_advance_mod.F90 5! Generated at: 2015-04-12 19:37:49 6! KGEN version: 0.4.9 7 8 9 10 MODULE prim_advance_mod 11 USE kgen_utils_mod, ONLY : kgen_dp, check_t, kgen_init_check, kgen_print_check 12 USE element_mod, ONLY : kgen_read_mod9 => kgen_read 13 USE element_mod, ONLY : kgen_verify_mod9 => kgen_verify 14 ! _EXTERNAL 15 IMPLICIT NONE 16 PRIVATE 17 PUBLIC compute_and_apply_rhs 18 CONTAINS 19 20 ! write subroutines 21 ! No subroutines 22 ! No module extern variables 23 24 25 26 27 28 29 30 31 32 ! 33 ! phl notes: output is stored in first argument. Advances from 2nd argument using tendencies evaluated at 3rd rgument: 34 ! phl: for offline winds use time at 3rd argument (same as rhs currently) 35 ! 36 37 SUBROUTINE compute_and_apply_rhs(elem, kgen_unit) 38 USE kgen_utils_mod, ONLY : kgen_dp, check_t, kgen_init_check, kgen_print_check 39 ! =================================== 40 ! compute the RHS, accumulate into u(np1) and apply DSS 41 ! 42 ! u(np1) = u(nm1) + dt2*DSS[ RHS(u(n0)) ] 43 ! 44 ! This subroutine is normally called to compute a leapfrog timestep 45 ! but by adjusting np1,nm1,n0 and dt2, many other timesteps can be 46 ! accomodated. For example, setting nm1=np1=n0 this routine will 47 ! take a forward euler step, overwriting the input with the output. 48 ! 49 ! qn0 = timelevel used to access Qdp() in order to compute virtual Temperature 50 ! qn0=-1 for the dry case 51 ! 52 ! if dt2<0, then the DSS'd RHS is returned in timelevel np1 53 ! 54 ! Combining the RHS and DSS pack operation in one routine 55 ! allows us to fuse these two loops for more cache reuse 56 ! 57 ! Combining the dt advance and DSS unpack operation in one routine 58 ! allows us to fuse these two loops for more cache reuse 59 ! 60 ! note: for prescribed velocity case, velocity will be computed at 61 ! "real_time", which should be the time of timelevel n0. 62 ! 63 ! 64 ! =================================== 65 USE kinds, ONLY: real_kind 66 USE dimensions_mod, ONLY: np 67 USE dimensions_mod, ONLY: nlev 68 USE element_mod, ONLY: element_t 69 USE prim_si_mod, ONLY: preq_hydrostatic 70 IMPLICIT NONE 71 integer, intent(in) :: kgen_unit 72 INTEGER*8 :: kgen_intvar, start_clock, stop_clock, rate_clock 73 TYPE(check_t):: check_status 74 REAL(KIND=kgen_dp) :: tolerance 75 TYPE(element_t), intent(inout), target :: elem(:) 76 ! weighting for eta_dot_dpdn mean flux 77 ! local 78 ! surface pressure for current tiime level 79 REAL(KIND=real_kind), pointer, dimension(:,:,:) :: phi 80 REAL(KIND=real_kind), pointer :: ref_phi(:,:,:) => NULL() 81 REAL(KIND=real_kind), dimension(np,np,nlev) :: t_v 82 ! half level vertical velocity on p-grid 83 ! temporary field 84 ! generic gradient storage 85 ! v.grad(T) 86 ! kinetic energy + PHI term 87 ! lat-lon coord version 88 ! vorticity 89 REAL(KIND=real_kind), dimension(np,np,nlev) :: p ! pressure 90 REAL(KIND=real_kind), dimension(np,np,nlev) :: dp ! delta pressure 91 ! inverse of delta pressure 92 ! temperature vertical advection 93 ! v.grad(p) 94 ! half level pressures on p-grid 95 ! velocity vertical advection 96 INTEGER :: ie 97 !JMD call t_barrierf('sync_compute_and_apply_rhs', hybrid%par%comm) 98 tolerance = 1.E-14 99 CALL kgen_init_check(check_status, tolerance) 100 CALL kgen_read_real_real_kind_dim3_ptr(phi, kgen_unit) 101 READ(UNIT=kgen_unit) t_v 102 READ(UNIT=kgen_unit) p 103 READ(UNIT=kgen_unit) dp 104 READ(UNIT=kgen_unit) ie 105 106 CALL kgen_read_real_real_kind_dim3_ptr(ref_phi, kgen_unit) 107 108 109 ! call to kernel 110 CALL preq_hydrostatic(phi, elem(ie)%state%phis, t_v, p, dp) 111 ! kernel verification for output variables 112 CALL kgen_verify_real_real_kind_dim3_ptr( "phi", check_status, phi, ref_phi) 113 CALL kgen_print_check("preq_hydrostatic", check_status) 114 CALL system_clock(start_clock, rate_clock) 115 DO kgen_intvar=1,10 116 CALL preq_hydrostatic(phi, elem(ie) % state % phis, t_v, p, dp) 117 END DO 118 CALL system_clock(stop_clock, rate_clock) 119 WRITE(*,*) 120 PRINT *, "Elapsed time (sec): ", (stop_clock - start_clock)/REAL(rate_clock*10) 121 ! ============================================================= 122 ! Insert communications here: for shared memory, just a single 123 ! sync is required 124 ! ============================================================= 125 CONTAINS 126 127 ! write subroutines 128 SUBROUTINE kgen_read_real_real_kind_dim3_ptr(var, kgen_unit, printvar) 129 INTEGER, INTENT(IN) :: kgen_unit 130 CHARACTER(*), INTENT(IN), OPTIONAL :: printvar 131 real(KIND=real_kind), INTENT(OUT), POINTER, DIMENSION(:,:,:) :: var 132 LOGICAL :: is_true 133 INTEGER :: idx1,idx2,idx3 134 INTEGER, DIMENSION(2,3) :: kgen_bound 135 136 READ(UNIT = kgen_unit) is_true 137 138 IF ( is_true ) THEN 139 READ(UNIT = kgen_unit) kgen_bound(1, 1) 140 READ(UNIT = kgen_unit) kgen_bound(2, 1) 141 READ(UNIT = kgen_unit) kgen_bound(1, 2) 142 READ(UNIT = kgen_unit) kgen_bound(2, 2) 143 READ(UNIT = kgen_unit) kgen_bound(1, 3) 144 READ(UNIT = kgen_unit) kgen_bound(2, 3) 145 ALLOCATE(var(kgen_bound(2, 1) - kgen_bound(1, 1) + 1, kgen_bound(2, 2) - kgen_bound(1, 2) + 1, kgen_bound(2, 3) - kgen_bound(1, 3) + 1)) 146 READ(UNIT = kgen_unit) var 147 IF ( PRESENT(printvar) ) THEN 148 PRINT *, "** " // printvar // " **", var 149 END IF 150 END IF 151 END SUBROUTINE kgen_read_real_real_kind_dim3_ptr 152 153 154 ! verify subroutines 155 SUBROUTINE kgen_verify_real_real_kind_dim3_ptr( varname, check_status, var, ref_var) 156 character(*), intent(in) :: varname 157 type(check_t), intent(inout) :: check_status 158 real(KIND=real_kind), intent(in), DIMENSION(:,:,:), POINTER :: var, ref_var 159 real(KIND=real_kind) :: nrmsdiff, rmsdiff 160 real(KIND=real_kind), allocatable, DIMENSION(:,:,:) :: temp, temp2 161 integer :: n 162 IF ( ASSOCIATED(var) ) THEN 163 check_status%numTotal = check_status%numTotal + 1 164 IF ( ALL( var == ref_var ) ) THEN 165 166 check_status%numIdentical = check_status%numIdentical + 1 167 if(check_status%verboseLevel > 1) then 168 WRITE(*,*) 169 WRITE(*,*) "All elements of ", trim(adjustl(varname)), " are IDENTICAL." 170 !WRITE(*,*) "KERNEL: ", var 171 !WRITE(*,*) "REF. : ", ref_var 172 IF ( ALL( var == 0 ) ) THEN 173 if(check_status%verboseLevel > 2) then 174 WRITE(*,*) "All values are zero." 175 end if 176 END IF 177 end if 178 ELSE 179 allocate(temp(SIZE(var,dim=1),SIZE(var,dim=2),SIZE(var,dim=3))) 180 allocate(temp2(SIZE(var,dim=1),SIZE(var,dim=2),SIZE(var,dim=3))) 181 182 n = count(var/=ref_var) 183 where(abs(ref_var) > check_status%minvalue) 184 temp = ((var-ref_var)/ref_var)**2 185 temp2 = (var-ref_var)**2 186 elsewhere 187 temp = (var-ref_var)**2 188 temp2 = temp 189 endwhere 190 nrmsdiff = sqrt(sum(temp)/real(n)) 191 rmsdiff = sqrt(sum(temp2)/real(n)) 192 193 if(check_status%verboseLevel > 0) then 194 WRITE(*,*) 195 WRITE(*,*) trim(adjustl(varname)), " is NOT IDENTICAL." 196 WRITE(*,*) count( var /= ref_var), " of ", size( var ), " elements are different." 197 if(check_status%verboseLevel > 1) then 198 WRITE(*,*) "Average - kernel ", sum(var)/real(size(var)) 199 WRITE(*,*) "Average - reference ", sum(ref_var)/real(size(ref_var)) 200 endif 201 WRITE(*,*) "RMS of difference is ",rmsdiff 202 WRITE(*,*) "Normalized RMS of difference is ",nrmsdiff 203 end if 204 205 if (nrmsdiff > check_status%tolerance) then 206 check_status%numFatal = check_status%numFatal+1 207 else 208 check_status%numWarning = check_status%numWarning+1 209 endif 210 211 deallocate(temp,temp2) 212 END IF 213 END IF 214 END SUBROUTINE kgen_verify_real_real_kind_dim3_ptr 215 216 END SUBROUTINE compute_and_apply_rhs 217 !TRILINOS 218 219 220 END MODULE prim_advance_mod 221