ROMS
Loading...
Searching...
No Matches
oyster_floats_inp.h
Go to the documentation of this file.
1 SUBROUTINE read_fltbiopar (model, inp, out, Lwrite)
2!
3!git $Id$
4!================================================== Hernan G. Arango ===
5! Copyright (c) 2002-2025 The ROMS Group !
6! Licensed under a MIT/X style license !
7! See License_ROMS.md !
8!=======================================================================
9! !
10! This routine reads in input biological floats parameters. !
11! !
12!=======================================================================
13!
14 USE mod_param
15 USE mod_parallel
16 USE mod_behavior
17 USE mod_iounits
18 USE mod_ncparam
19 USE mod_scalars
20!
22!
23 implicit none
24!
25! Imported variable declarations
26!
27 logical, intent(in) :: Lwrite
28 integer, intent(in) :: model, inp, out
29!
30! Local variable declarations.
31!
32 integer :: Npts, Nval
33 integer :: i, j, igrid, mc, nc, ng, status
34 integer :: Ivalue(1)
35
36 real(r8) :: Rvalue(1)
37
38 real(dp), dimension(nRval) :: Rval
39
40 character (len=35) :: frmt
41 character (len=40) :: KeyWord
42 character (len=256) :: line
43 character (len=256), dimension(nCval) :: Cval
44
45 character (len=1 ), parameter :: blank = ' '
46!
47!-----------------------------------------------------------------------
48! Read in initial float locations.
49!-----------------------------------------------------------------------
50!
51! Allocate oyster model parameter.
52!
54!
55! Notice I added one when allocating local scratch arrays to avoid
56! out of bounds in some compilers when reading the last blank line
57! which signal termination of input data.
58!
59 DO WHILE (.true.)
60 READ (inp,'(a)',err=10,END=30) line
61 status=decode_line(line, keyword, nval, cval, rval)
62 IF (status.gt.0) THEN
63 SELECT CASE (trim(keyword))
64 CASE ('Larvae_size0')
65 npts=load_r(nval, rval, ngrids, larvae_size0)
66 CASE ('Larvae_GR0')
67 npts=load_r(nval, rval, ngrids, larvae_gr0)
68 CASE ('settle_size')
69 npts=load_r(nval, rval, ngrids, settle_size)
70 CASE ('food_supply')
71 npts=load_r(nval, rval, ngrids, food_supply)
72 CASE ('turb_ambi')
73 npts=load_r(nval, rval, ngrids, turb_ambi)
74 CASE ('turb_crit')
75 npts=load_r(nval, rval, ngrids, turb_crit)
76 CASE ('turb_slop')
77 npts=load_r(nval, rval, ngrids, turb_slop)
78 CASE ('turb_axis')
79 npts=load_r(nval, rval, ngrids, turb_axis)
80 CASE ('turb_base')
81 npts=load_r(nval, rval, ngrids, turb_base)
82 CASE ('turb_rate')
83 npts=load_r(nval, rval, ngrids, turb_rate)
84 CASE ('turb_mean')
85 npts=load_r(nval, rval, ngrids, turb_mean)
86 CASE ('turb_size')
87 npts=load_r(nval, rval, ngrids, turb_size)
88 CASE ('swim_Tmin')
89 npts=load_r(nval, rval, ngrids, swim_tmin)
90 CASE ('swim_Tmax')
91 npts=load_r(nval, rval, ngrids, swim_tmax)
92 CASE ('swim_Sinc')
93 npts=load_r(nval, rval, ngrids, swim_sinc)
94 CASE ('swim_Sdec')
95 npts=load_r(nval, rval, ngrids, swim_sdec)
96 CASE ('slope_Sinc')
97 npts=load_r(nval, rval, ngrids, slope_sinc)
98 CASE ('slope_Sdec')
99 npts=load_r(nval, rval, ngrids, slope_sdec)
100 CASE ('sink_base')
101 npts=load_r(nval, rval, ngrids, sink_base)
102 CASE ('sink_rate')
103 npts=load_r(nval, rval, ngrids, sink_rate)
104 CASE ('sink_size')
105 npts=load_r(nval, rval, ngrids, sink_size)
106 CASE ('swim_Im')
107 npts=load_i(nval, rval, 1, ivalue)
108 swim_im=ivalue(1)
109 CASE ('swim_Jm')
110 npts=load_i(nval, rval, 1, ivalue)
111 swim_jm=ivalue(1)
112 CASE ('swim_L0')
113 npts=load_r(nval, rval, 1, rvalue)
114 swim_l0=rvalue(1)
115 CASE ('swim_T0')
116 npts=load_r(nval, rval, 1, rvalue)
117 swim_t0=rvalue(1)
118 CASE ('swim_DL')
119 npts=load_r(nval, rval, 1, rvalue)
120 swim_dl=rvalue(1)
121 CASE ('swim_DT')
122 npts=load_r(nval, rval, 1, rvalue)
123 swim_dt=rvalue(1)
124 CASE ('swim_table')
125 IF (.not.allocated(swim_table)) THEN
126 allocate ( swim_table(swim_im,swim_jm) )
127 swim_table=0.0_r8
128 dmem(1)=dmem(1)+real(swim_im*swim_jm,r8)
129 END IF
130 READ (inp,*,err=20,END=30) &
131 ((swim_table(i,j),i=1,swim_im),j=1,swim_jm)
132 CASE ('Gfactor_Im')
133 npts=load_i(nval, rval, 1, ivalue)
134 gfactor_im=ivalue(1)
135 CASE ('Gfactor_Jm')
136 npts=load_i(nval, rval, 1, ivalue)
137 gfactor_jm=ivalue(1)
138 CASE ('Gfactor_S0')
139 npts=load_r(nval, rval, 1, rvalue)
140 gfactor_s0=rvalue(1)
141 CASE ('Gfactor_T0')
142 npts=load_r(nval, rval, 1, rvalue)
143 gfactor_t0=rvalue(1)
144 CASE ('Gfactor_DS')
145 npts=load_r(nval, rval, 1, rvalue)
146 gfactor_ds=rvalue(1)
147 CASE ('Gfactor_DT')
148 npts=load_r(nval, rval, 1, rvalue)
149 gfactor_dt=rvalue(1)
150 CASE ('Gfactor_table')
151 IF (.not.allocated(gfactor_table)) THEN
152 allocate ( gfactor_table(gfactor_im,gfactor_jm) )
153 gfactor_table=0.0_r8
154 dmem(1)=dmem(1)+real(gfactor_im*gfactor_jm)
155 END IF
156 READ (inp,*,err=20,END=30) &
157 ((gfactor_table(i,j),i=1,gfactor_im),j=1,gfactor_jm)
158 CASE ('Grate_Im')
159 npts=load_i(nval, rval, 1, ivalue)
160 grate_im=ivalue(1)
161 CASE ('Grate_Jm')
162 npts=load_i(nval, rval, 1, ivalue)
163 grate_jm=ivalue(1)
164 CASE ('Grate_F0')
165 npts=load_r(nval, rval, 1, rvalue)
166 grate_f0=rvalue(1)
167 CASE ('Grate_L0')
168 npts=load_r(nval, rval, 1, rvalue)
169 grate_l0=rvalue(1)
170 CASE ('Grate_DF')
171 npts=load_r(nval, rval, 1, rvalue)
172 grate_df=rvalue(1)
173 CASE ('Grate_DL')
174 npts=load_r(nval, rval, 1, rvalue)
175 grate_dl=rvalue(1)
176 CASE ('Grate_table')
177 IF (.not.allocated(grate_table)) THEN
178 allocate ( grate_table(grate_im,grate_jm) )
179 grate_table=0.0_r8
180 dmem(1)=dmem(1)+real(grate_im*grate_jm,r8)
181 END IF
182 READ (inp,*,err=20,END=30) &
183 ((grate_table(i,j),i=1,grate_im),j=1,grate_jm)
184 END SELECT
185 END IF
186 END DO
187 10 IF (master) WRITE (out,40) line
188 exit_flag=4
189 RETURN
190 20 IF (master) WRITE (out,50) trim(keyword)
191 exit_flag=4
192 RETURN
193 30 CLOSE (inp)
194!
195!-----------------------------------------------------------------------
196! Report input parameters.
197!-----------------------------------------------------------------------
198!
199 IF (master.and.lwrite) THEN
200 DO ng=1,ngrids
201 WRITE (out,60) ng
202 WRITE (out,70) larvae_size0(ng), 'Larvae_size0', &
203 & 'Initial larvae size (um).'
204 WRITE (out,70) larvae_gr0(ng), 'Larvae_GR0', &
205 & 'Initial larvae growth rate (um/day).'
206 WRITE (out,70) settle_size(ng), 'settle_size', &
207 & 'Larvae settlement size (um).'
208 WRITE (out,70) food_supply(ng), 'food_supply', &
209 & 'Food supply (mg Carbon/l).'
210 WRITE (out,70) turb_ambi(ng), 'turb_ambi', &
211 & 'Ambient turbidity level (g/l).'
212 WRITE (out,70) turb_crit(ng), 'turb_crit', &
213 & 'Critical turbidity value (g/l).'
214 WRITE (out,70) turb_slop(ng), 'turb_slop', &
215 & 'Turbidity linear slope (l/g).'
216 WRITE (out,70) turb_axis(ng), 'turb_axis', &
217 & 'Turbidity linear axis crossing (g/l).'
218 WRITE (out,70) turb_base(ng), 'turb_base', &
219 & 'Turbidity exponential base factor (g/l).'
220 WRITE (out,70) turb_rate(ng), 'turb_rate', &
221 & 'Turbidity exponential rate (l/g).'
222 WRITE (out,70) turb_mean(ng), 'turb_mean', &
223 & 'Turbidity exponential mean (g/l).'
224 WRITE (out,70) turb_size(ng), 'turb_size', &
225 & 'Minimum larvae size (um) affected by turbidity.'
226 WRITE (out,70) swim_tmin(ng), 'swim_Tmin', &
227 & 'Minimum swimming time fraction.'
228 WRITE (out,70) swim_tmax(ng), 'swim_Tmax', &
229 & 'Maximum swimming time fraction.'
230 WRITE (out,70) swim_sinc(ng), 'swim_Sinc', &
231 & 'Swimming, active fraction due to increasing salinity.'
232 WRITE (out,70) swim_sdec(ng), 'swim_Sdec', &
233 & 'Swimming, active fraction due to decreasing salinity.'
234 WRITE (out,70) slope_sinc(ng), 'slope_Sinc', &
235 & 'Swimming, coefficient due to increasing salinity.'
236 WRITE (out,70) slope_sdec(ng), 'slope_Sdec', &
237 & 'Swimming, coefficient due to increasing salinity.'
238 WRITE (out,70) sink_base(ng), 'sink_base', &
239 & 'Sinking, exponential base factor (mm/s).'
240 WRITE (out,70) sink_rate(ng), 'sink_rate', &
241 & 'Sinking, exponential rate factor (1/um).'
242 WRITE (out,70) sink_size(ng), 'sink_mean', &
243 & 'Sinking, exponential mean size (um).'
244 WRITE (out,80) swim_im, 'swim_Im', &
245 & 'Swim table, number of values in larval size I-axis.'
246 WRITE (out,80) swim_jm, 'swim_Jm', &
247 & 'Swim table, number of values in temperature J-axis.'
248 WRITE (out,70) swim_l0, 'swim_L0', &
249 & 'Swim table, starting value for larval size I-axis.'
250 WRITE (out,70) swim_t0, 'swim_T0', &
251 & 'Swim table, starting value for temperature J-axis.'
252 WRITE (out,70) swim_dl, 'swim_DL', &
253 & 'Swim table, larval size I-axis increment.'
254 WRITE (out,70) swim_dt, 'swim_DT', &
255 & 'Swim table, temperature J-axis increment.'
256 WRITE (out,80) gfactor_im, 'Gfactor_Im', &
257 & 'Gfactor table, number of values in salinity I-axis.'
258 WRITE (out,80) gfactor_jm, 'Gfactor_Jm', &
259 & 'Gfactor table, number of values in temperature J-axis.'
260 WRITE (out,70) gfactor_s0, 'Gfactor_S0', &
261 & 'Gfactor table, starting value for salinity I-axis.'
262 WRITE (out,70) gfactor_t0, 'Gfactor_T0', &
263 & 'Gfactor table, starting value for temperature J-axis.'
264 WRITE (out,70) gfactor_ds, 'Gfactor_DS', &
265 & 'Gfactor table, starting value for salinity I-axis.'
266 WRITE (out,70) gfactor_dt, 'Gfactor_DT', &
267 & 'Gfactor table, starting value for temperature J-axis.'
268 WRITE (out,80) grate_im, 'Grate_Im', &
269 & 'Grate table, number of values in food supply I-axis.'
270 WRITE (out,80) grate_jm, 'Grate_Jm', &
271 & 'Grate table, number of values in larval size J-axis.'
272 WRITE (out,70) grate_f0, 'Grate_F0', &
273 & 'Grate table, starting value for food supply I-axis.'
274 WRITE (out,70) grate_l0, 'Grate_L0', &
275 & 'Grate table, starting value for larval size J-axis.'
276 WRITE (out,70) grate_df, 'Grate_DF', &
277 & 'Grate table, food supply I-axis increment.'
278 WRITE (out,70) grate_dl, 'Grate_DL', &
279 & 'Grate table, larval size J-axis increment.'
280 END DO
281 END IF
282
283 40 FORMAT (/,' READ_FloatsBioPar - Error while processing line: ',/, &
284 & a)
285 50 FORMAT (/,' READ_FloatsBioPar - Error reading look table: ',a)
286 60 FORMAT (/,/,' Biological Floats Behavior Parameters, Grid: ',i2.2,&
287 & /, ' ===============================================',/)
288 70 FORMAT (1p,e11.4,2x,a,t32,a)
289 80 FORMAT (1x,i10,2x,a,t32,a)
290
291 RETURN
292 END SUBROUTINE read_fltbiopar
subroutine allocate_behavior
type(t_io), dimension(:), allocatable err
subroutine read_fltbiopar(model, inp, out, lwrite)