source: trunk/LMDZ.GENERIC/libf/phystd/check_fields.F90 @ 2987

Last change on this file since 2987 was 2873, checked in by jleconte, 23 months ago

Correction in check fields

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