source: trunk/LMDZ.MARS/libf/phymars/check_fields.F90 @ 2800

Last change on this file since 2800 was 2570, checked in by emillour, 3 years ago

Mars GCM:
Improve tests by "check_physics_fields": also look for
negative values of tracers.
EM

File size: 6.5 KB
Line 
1module check_fields_mod
2
3  real,parameter :: default_temp_min=50.  ! minimum reasonable temperature (K)
4  real,parameter :: default_temp_max=350. ! maximum reasonable temperature (K)
5
6  real,parameter :: default_wind_max=500. ! maximum reasonable wind magnitude (m/s)
7
8  real,parameter :: default_ps_min=80.  ! minimum reasonable surface pressure (Pa)
9  real,parameter :: default_ps_max=2000. ! maximum reasonable surface pressure (Pa)
10
11contains
12
13subroutine check_physics_fields(message,temp,u,v,pplev,q)
14use dimphy, only: klon, klev
15use tracer_mod, only: nqmx
16implicit none
17character(len=*),intent(in):: message ! text message for outputs
18real,intent(in) :: temp(klon,klev)
19real,intent(in) :: u(klon,klev) ! zonal wind (m/s)
20real,intent(in) :: v(klon,klev) ! meridional wind (m/s)
21real,intent(in) :: pplev(klon,klev+1) ! pressure at level interfaces (Pa)
22real,intent(in) :: q(klon,klev,nqmx) ! tracer mixing ratio (.../kg_air)
23
24character(len=50) :: name="check_physics_fields"
25logical :: ok_t,ok_w,ok_ps,ok_q
26
27! 1. Initialisations
28ok_q=.true.
29
30! 2. Check temperature, winds and surface pressure
31call check_temperature(message,temp,ok_t)
32call check_winds(message,u,v,ok_w)
33call check_ps(message,pplev(:,1),ok_ps)
34if (nqmx>=1) then
35  call check_tracers(message,q,ok_q)
36endif
37
38if ((.not.ok_t).or.(.not.ok_w).or.(.not.ok_ps).or.(.not.ok_q)) then
39  ! Something is wrong, might as well stop here
40  call abort_physic(trim(name),trim(message)//" Invalid field values",1)
41endif
42
43end subroutine check_physics_fields
44
45
46subroutine check_temperature(message,temp,ok,temp_min,temp_max)
47use dimphy, only: klon, klev
48implicit none
49character(len=*),intent(in):: message ! text message for outputs
50real,intent(in) :: temp(klon,klev)
51logical,intent(out) :: ok ! returns .true. if everything OK, .false. otherwise
52real,intent(in),optional :: temp_min ! user provided acceptable minimum
53real,intent(in),optional :: temp_max ! user provided acceptable maximum
54
55character(len=50) :: name="check_temperature"
56real :: tmin,tmax
57integer :: i,k
58
59! 0. Check optional inputs
60if (present(temp_min)) then
61  tmin=temp_min
62else
63  tmin=default_temp_min
64endif
65
66if (present(temp_max)) then
67  tmax=temp_max
68else
69  tmax=default_temp_max
70endif
71
72! 1. initializations
73ok=.true.
74
75! 2. do the checking
76do i=1,klon
77  do k=1,klev
78    ! look for NaN
79    if (temp(i,k).ne.temp(i,k)) then
80      ok=.false.
81      write(*,*)trim(message)//" "//trim(name)//" temp(i,k)=",temp(i,k),&
82                " for i=",i," k=",k
83    endif
84    ! check for temperatures too low
85    if (temp(i,k).lt.tmin) then
86      ok=.false.
87      write(*,*)trim(message)//" "//trim(name)//" temp(i,k)=",temp(i,k),&
88      " for i=",i," k=",k," <",tmin
89    endif
90    ! check for temperatures too high
91    if (temp(i,k).gt.tmax) then
92      ok=.false.
93      write(*,*)trim(message)//" "//trim(name)//" temp(i,k)=",temp(i,k),&
94      " for i=",i," k=",k," >",tmax
95    endif
96  enddo
97enddo
98
99end subroutine check_temperature
100
101subroutine check_winds(message,u,v,ok,wind_max)
102use dimphy, only: klon, klev
103implicit none
104character(len=*),intent(in):: message ! text message for outputs
105real,intent(in) :: u(klon,klev) ! zonal wind (m/s)
106real,intent(in) :: v(klon,klev) ! meridional wind (m/s)
107logical,intent(out) :: ok ! returns .true. if everything OK, .false. otherwise
108real,intent(in),optional :: wind_max ! user provided acceptable maximum magnitude
109
110character(len=50) :: name="check_winds"
111real :: wmax
112integer :: i,k
113
114! 0. Check optional inputs
115
116if (present(wind_max)) then
117  wmax=wind_max
118else
119  wmax=default_wind_max
120endif
121
122! 1. initializations
123ok=.true.
124
125! 2. do the checking
126do i=1,klon
127  do k=1,klev
128    ! look for NaN
129    if (u(i,k).ne.u(i,k)) then
130      ok=.false.
131      write(*,*)trim(message)//" "//trim(name)//" u(i,k)=",u(i,k),&
132      " for i=",i," k=",k
133    endif
134    if (v(i,k).ne.v(i,k)) then
135      ok=.false.
136      write(*,*)trim(message)//" "//trim(name)//" v(i,k)=",v(i,k),&
137      " for i=",i," k=",k
138    endif
139    ! check for magnitudes too high
140    if (abs(u(i,k)).gt.wmax) then
141      ok=.false.
142      write(*,*)trim(message)//" "//trim(name)//" u(i,k)=",u(i,k),&
143      " for i=",i," k=",k," >",wmax
144    endif
145    if (abs(v(i,k)).gt.wmax) then
146      ok=.false.
147      write(*,*)trim(message)//" "//trim(name)//" v(i,k)=",v(i,k),&
148      " for i=",i," k=",k," >",wmax
149    endif
150  enddo
151enddo
152
153end subroutine check_winds
154
155subroutine check_ps(message,ps,ok,ps_min,ps_max)
156use dimphy, only: klon
157implicit none
158character(len=*),intent(in):: message ! text message for outputs
159real,intent(in) :: ps(klon)
160logical,intent(out) :: ok ! returns .true. if everything OK, .false. otherwise
161real,intent(in),optional :: ps_min ! user provided acceptable minimum
162real,intent(in),optional :: ps_max ! user provided acceptable maximum
163
164character(len=50) :: name="check_ps"
165real :: pmin,pmax
166integer :: i
167
168! 0. Check optional inputs
169if (present(ps_min)) then
170  pmin=ps_min
171else
172  pmin=default_ps_min
173endif
174
175if (present(ps_max)) then
176  pmax=ps_max
177else
178  pmax=default_ps_max
179endif
180
181! 1. initializations
182ok=.true.
183
184! 2. do the checking
185do i=1,klon
186  ! look for NaN
187  if (ps(i).ne.ps(i)) then
188    ok=.false.
189    write(*,*)trim(message)//" "//trim(name)//" ps(i)=",ps(i)," for i=",i
190  endif
191  ! check for pressures too low
192  if (ps(i).lt.pmin) then
193    ok=.false.
194    write(*,*)trim(message)//" "//trim(name)//" ps(i)=",ps(i)," for i=",i,&
195    " <",pmin
196  endif
197  ! check for pressures too high
198  if (ps(i).gt.pmax) then
199    ok=.false.
200    write(*,*)trim(message)//" "//trim(name)//" ps(i)=",ps(i)," for i=",i,&
201    " >",pmax
202  endif
203enddo
204
205end subroutine check_ps
206
207subroutine check_tracers(message,q,ok)
208use dimphy, only: klon, klev
209use tracer_mod, only: nqmx,noms
210implicit none
211character(len=*),intent(in):: message ! text message for outputs
212real,intent(in) :: q(klon,klev,nqmx) ! tracer mixing ratio
213logical,intent(inout) :: ok ! set to .false. if something is wrong
214
215character(len=50) :: name="check_tracers"
216integer :: i,k,iq
217integer :: nb_bad_neg ! number of problematic points (negative values)
218integer :: nb_bad_nan ! number of problematic points (NaNs)
219
220! 1. initializations
221nb_bad_nan=0
222
223! 2. do the checking
224do iq=1,nqmx
225 nb_bad_neg=0 ! initialize
226 do i=1,klon
227  do k=1,klev
228    ! look for NaN
229    if (q(i,k,iq).ne.q(i,k,iq)) then
230      ok=.false.
231      nb_bad_nan=nb_bad_nan+1
232    endif
233    ! look for negative values
234    if (q(i,k,iq).lt.0) then
235!      ok=.false.
236      nb_bad_neg=nb_bad_neg+1
237    endif
238  enddo
239 enddo
240 if (nb_bad_neg>0) then
241   write(*,*)trim(message)//" "//trim(name)//" tracer "//trim(noms(iq))//&
242   " contains ",nb_bad_neg," negative values!"
243 endif
244enddo !of do iq=1,nqmx
245
246end subroutine check_tracers
247
248end module check_fields_mod
Note: See TracBrowser for help on using the repository browser.