*FB2autoVarycbD.f Built on FB2auto.f - 9 Mar 00 * * Solves RGR(C-gain-limited) = RGR(C-use-limited) (the new PS inhibition * phenomena), to get RGR(C-limited) and RGR(C-limited) = RGR(N-limited), * The equations remain analytic but the solution is not in closed form, * since the solution of 1st equality gives fN appearing in square roots * in both numerator and denominator of RGR(C-limited), so that solving * 2nd FB equation for fN (thus, for RGR(FB) ) is not possible. * Alternative, simple method used here: do binary search in fN: fast, * simple, short program that any one can use and in which anyone can * see implications of changes in physiol. params., CO2 level, MR infection, * and soil diffusive limitations (in later version) * * Features over and above those in FB2.f: allows solution for RGR and fN * over a preset range of any one parameter, with corresponding calc. * of sensitivity param. and storage of RGR and fN results in a file * * Features in addition to those in FB2auto.f: allows range of either * cb or D in computing Vmax, and computes Sens(RGR, to Vmax) at each pt. * parameter (nparam=13) implicit real*8 (a-h,o-z) real*8 kC, Ke real*8 Km, Ja character*6 aname(nparam),cmd character*24 outfile,opened(50) real*8 value(nparam),valueprev(nparam) character*55 message1(nparam),message2(nparam)*10 character*55 definition(nparam),def2(nparam) dimension value0(nparam) equivalence (value(1),alphaL),(value(2),pstar),(value(3),r), + (value(4),beta),(value(5),vplant),(value(6),alphaM), + (value(7),vM),(value(8),kC),(value(9),Ke),(value(10),G), + (value(11),RGRmax),(value(12),fNlo0),(value(13),fNhi0) common CalphaL,Cpstar,Cr,Cbeta,Cvplant, CalphaM,CvM,CkC,CKe, + CG,CRGRmax,CfC data aname/'alphaL','pstar','r','beta','vplant', + 'alphaM','vM','kC','Ke','G','RGRmax','fNlo0','fNhi0'/ data message1/'none (leaf dry mass)/(shoot dry mass)', + 'g(PSate) gN^-1 d^-1 (i.e., per photoperiod)', + 'none (root dry mass)/(shoot dry mass)', + 'gDM/g(PSate', + 'gN g(DMroot)^-1 d^-1 (N uptake per unit root dry mass)', + 'g(MRfungi)/g(DMroot only)','gN g(MRfungi)^-1 d^-1', + 2*'molar conc. sucrose in tissue', + 'none (metab. rate of MRfungi per mass)/(same for plant)', + 'd^-1 (same as gDM gDM^-1 d^-1)','none (or gN gDM^-1)', + 'none (or gN gDM^-1)'/ data message2/'0.5','26.','0.4','0.5','0.017', + '0.3','0.025','1.','0.5','3.','0.3','0.01','0.06'/ data value0/0.5,26.,0.4,0.5,0.017,0.3,0.025,1.,0.5,3., + 0.3,0.01,0.06/ data definition/'frac. of shoot mass as leaf', + 'photosyn. N-use efficiency', + 'root:shoot ratio=(root mass)/(shoot mass)', + 'frac. of raw PSate converted to final biomass', + 'mean N-uptake rate per dry mass of root', + '(dry mass of MR fungi)/(dry mass of root only)', + 'mean N-uptake rate per dry mass of MR fungi', + 'tissue sucrose conc. for 1/2 inhib. of PS', + 'tissue sucrose conc. for 1/2 sat. of RGR', + '(metab. rate per mass of MR fungi)/(rate for plant)', + 'maximum RGR (sucrose-saturated)', + 'lower limit of fN in search for solution', + 'upper limit of fN in search for solution'/ data def2/3*'','(acctg. for biosyn. and maint. losses)',9*''/ * Keep track of results files already opened iopen=0 * Read in parameters for first run 60 write(6,'(/,"We will enter values for:",/)') do i=1,nparam write(6,'(a6," = ",a55)')aname(i),definition(i) if(aname(i).eq.'beta')write(6,'(10x,a55)')def2(i) enddo write(6,'(/,"Units will be noted ", + "as you are prompted for numerical values.",/,1x,"For vplant", + " and pstar, several options are given if you need to convert", + /,5x,"from different units in your own data",/)') do 500 i=1,nparam write(6,'("Enter value of ",a6," = ",a55,/,5x,"Units: ", + a55,/,5x,"Typical value: ",a10,3x,"Current value: ",g10.5)') + aname(i),definition(i),message1(i),message2(i),value0(i) if(aname(i).eq.'pstar')then call pstarread(pstar) elseif(aname(i).eq.'vplant')then call vplantread(vplant) elseif(aname(i).eq.'vM')then call vMread(vplant,vM) else read(5,*)value(i) ! If not a special case, just read in value endif 500 enddo ! End of do-loop for entering values of all parameters initially tol=0.0000005 en=dlog((fNhi0-fNlo0)/tol)/0.693 nbin=ifix(en) * Also put equivalenced variables into new variables that can be passed * in common - it is dangerous to pass equivalenced vars. to common! CalphaL=alphaL Cpstar=pstar Cr=r Cbeta=beta Cvplant=vplant CalphaM=alphaM CvM=vM CkC=kC CKe=Ke CG=G CRGRmax=RGRmax call solver(fN,RGR,fNlo0,fNhi0,nbin) RGRlast=RGR fNlast=fN do i=1,nparam valueprev(i)=value(i) enddo * Now allow changes of parameter values 90 write(6,'(/,"Enter ''1'' to change a single parameter value,", + /,3x,"enter ''2'' to enter a range to explore in any ", + "parameter")') read(5,*,iostat=ierr)nchoice if(ierr.ne.0)then write(6,'("Error (probably letter instead of number); ", + "try again")') goto 90 endif if(nchoice.eq.1)then nchanged=0 write(6,'("To modify a parameter value, enter its name (up to", + " 6 char., left-justified)",/,5x,"You can return to a param.", + " as often as you like, and revise",/,8x,"any number of them,", + " in any order.",/,5x,"An entry of ''run'' ", + " requests a run with current",/,8x,"param. values.",/,5x, + "An entry of ''see'' causes printout of all params. and their", + /,8x," current values.",/,5x,"An entry of ''all'' is a ",/,8x, + "request to", + "re-enter all param. values.",/,5x,"An entry of ''stop'' stops", + " the program.")') 87 write(6,'("Enter the command of choice (name, ''run'', ''see'',", + "''all'', or ''stop'')")') 871 read(5,*,iostat=ierr)cmd if(ierr.ne.0)then write(6,'("Error (probably number instead of letter);", + " try again")') goto 87 endif if(cmd.eq.'stop')then stop elseif(cmd.eq.'run')then goto 200 ! Run the simulation elseif(cmd.eq.'all')then ncha=0 goto 60 elseif(cmd.eq.'see')then do i=1,nparam write(6,'(a6," =",g10.5,2x,a55)')aname(i),value(i), + message1(i) enddo goto 87 else ! Test to see if a param. name is recognizable do i=1,nparam if(aname(i).eq.cmd)then write(6,'("Enter value of ",a6," = ",a55,/,5x,"Units: ", + a55,/,5x,"Typical value: ",a10,3x,"Current value: ", + g10.5)')aname(i),definition(i),message1(i),message2(i), + value(i) if(i.le.nparam-2)then nchanged=nchanged+1 ichanged=i * Increment the counter for number of changes in physiol. or soil params., * and the identity of the parameter (as its number, i) but exclude * changes in the final two params. that only define the search interval * in fN, as (fNlo0, fNhi0) endif if(aname(i).eq.'pstar')then call pstarread(pstar) anamechanged='pstar' valuenew=pstar goto 87 elseif(aname(i).eq.'vplant')then vplantold=vplant call vplantread(vplant) anamechanged='vplant' valuenew=vplant write(6,'("Does vM stay in proportion to vplant (1) or", + " does vM itself stay constant (2) ?")') 993 read(5,*,iostat=ierr)iv if(ierr.ne.0)then write(6,'("Error (probably text instead of number);", + " try again")') goto 993 endif if(iv.eq.1)then vM=vM*vplant/vplantold endif goto 87 endif read(5,*)value(i) anamechanged=aname(i) valuenew=value(i) goto 87 endif enddo write(6,'("Command or param. name not recognized; try again", + )') goto 87 endif 200 tol=0.0000005 en=dlog((fNhi0-fNlo0)/tol)/0.693 nbin=ifix(en) * Also put equivalenced variables into new variables that can be passed * in common - it is dangerous to pass equivalenced vars. to common! CalphaL=alphaL Cpstar=pstar Cr=r Cbeta=beta Cvplant=vplant CalphaM=alphaM CvM=vM CkC=kC CKe=Ke CG=G CRGRmax=RGRmax call solver(fN,RGR,fNlo0,fNhi0,nbin) dRGR=RGR-RGRlast write(6,'("Change in RGR from last run is ",f8.5)')dRGR if((nchanged.gt.0).and.(nchanged.lt.2))then SensRGR=valueprev(ichanged)*dRGR/(RGRlast* + (value(ichanged)-valueprev(ichanged))) write(6,'("Only one plant or soil param. was changed.",/, + "Thus, we can calculate sensitivity of RGR to this", + " factor, defining it as",/,5x,"(relative change in RGR)/", + "(relative change in factor)",/,4x"(that is, d(ln RGR)/", + "d(ln factor). Example: a 10% change in vplant",/,4x, + "causes a 3.2% change in RGR (say, from 0.2 to 0.2064).", + " This is",/,4x,"0.32 of 10%, so the sensitivity is 0.32.", + /,1x,"The change in the factor ",a6," was from ",g10.5, + " to ",g10.5,/,5x,"so that Sens.=",f10.5)')aname(ichanged), + valueprev(ichanged),value(ichanged),SensRGR * Do same for sensitivity of fN to change in parameter value SensfN=valueprev(ichanged)*(fN-fNlast)/(fNlast* + (value(ichanged)-valueprev(ichanged))) write(6,'("Similarly, fN changed from ",g10.5," to ", + g10.5,/,5x," as a result of the change in ",a6,". The", + " sensitivity factor for fN,",/,5x,"defined similarly to", + " sensitivity of RGR, is ",f10.5)')fNlast,fN,aname(ichanged), + SensfN else write(6,'("Zero, or two or more parameters were changed;", + " we can''t calculate a",/,5x,"sensitivity of RGR to one", + " param.")') endif RGRlast=RGR fNlast=fN * Store current values in record of previous values, for next iteration; * the reason we update ALL param. values, and not just the one recorded * as being changed most recently, is that 2 or more parameters may have * been changed do i=1,nparam valueprev(i)=value(i) enddo goto 90 elseif(nchoice.eq.2)then itry=0 ! Haven't yet tried to open a file 101 write(6,'("Enter name of parameter to be varied over a range")') read(5,*)cmd 105 if(itry.eq.1)then write(6,'("Error; perhaps file exists already, but the", + "program has just started",/,5x,"and has no memory", + " of this; delete or move file and try again")') itry=0 ! Reset this for next round endif write(6,'("Enter name of file (up to 24 char. long) in which", + " to store results as ",/,10x, + "value,RGR,fN,SensRGR,SensfN;",/,5x,"If name is exactly ", + "''blank'', no file will be written",/,5x,"File will be in", + " ''old'' (append) mode (won''t overwrite prev. results", + /,8x,"in file of same name created earlier")') read(5,*)outfile itry=1 ! Got past first attempt to open file if(outfile.eq.'blank')then nowrite=1 else nowrite=0 * Test if this is a previously opened file if(iopen.gt.0)then do i=1,iopen if(outfile.eq.opened(i))then open(8,file=outfile,status='old') goto 110 endif enddo * Not found; must be a new file iopen=iopen+1 opened(iopen)=outfile open(8,file=outfile,status='new',err=105) else iopen=1 opened(1)=outfile open(8,file=outfile,status='new',err=105) endif ! End test for newness or oldness of file endif ! End test for any file to be opened (vs 'blank') * See if it's a recognizable name 110 if(cmd.eq.'cbandD')then i=5 ! identity of vplant as a parameter write(6,'("Enter a(m), b(m), Km(mol m-3), Vmaxa(mol ", + "m-2 s-1), fDMroot, Pper(h)")') read(5,*)a,b,Km,Vmaxa,fDMroot,Pper write(6,'("Enter D0, Df, dD; cb0, cbf, dcb")') read(5,*)D0,Df,dD,cb0,cbf,dcb do D=D0,Df,dD x=Vmaxa*a*dlog(b/a)/D do cb=cb0,cbf,dcb bbig=x-cb+Km ca=(-bbig+dsqrt(bbig*bbig+4.*Km*cb))*0.5 Ja=Vmaxa*ca/(ca+Km) ! This is mol m-2 s-1 fac=2.*ca+bbig dVm=1.-Km*x/((ca+Km)*fac) Cvplant=Ja*14.*(2./a)*(1.e-6/fDMroot)*Pper*3600. call solver(fN,RGR,fNlo0,fNhi0,nbin) value(i)=Cvplant SensRGR=valueprev(i)*(RGR-RGRlast)/(RGRlast* + (value(i)-valueprev(i))) SensRGRVmax=SensRGR*dVm valueprev(i)=Cvplant RGRlast=RGR * if(nowrite.eq.0)write(8,'(f10.5,10Pf10.5,0Pf10.5,3f10.5)') * + cb,D,RGR,SensRGR,dVm,SensRGRVmax if(nowrite.eq.0)write(8,'(f10.5,10Pf10.5,0Pf10.5,f10.5)') + cb,D,RGR,SensRGRVmax enddo enddo goto 90 endif do i=1,nparam if(aname(i).eq.cmd)then 111 write(6,'("Units of ",a6," must be standard units,",/,5x + a55,/,"Enter low value, high value, and step")')cmd, + message1(i) read(5,*,iostat=ierr)vlo,vhi,dval if(ierr.ne.0)then write(6,'("Error; try again (use spaces or commas as ", + "separators; use valid numbers)")') goto 111 endif do val=vlo,vhi,dval value(i)=val CalphaL=alphaL Cpstar=pstar Cr=r Cbeta=beta Cvplant=vplant CalphaM=alphaM CvM=vM CkC=kC CKe=Ke CG=G CRGRmax=RGRmax call solver(fN,RGR,fNlo0,fNhi0,nbin) SensRGR=valueprev(i)*(RGR-RGRlast)/(RGRlast* + (value(i)-valueprev(i))) SensfN=valueprev(i)*(fN-fNlast)/(fNlast* + (value(i)-valueprev(i))) write(6,'(a6,"= ",f10.5,"; Sens. of RGR to ",a6, + " is ",f7.4,"; Sens. of fN is ",f7.4)')aname(i), + val,aname(i),SensRGR,SensfN valueprev(i)=val RGRlast=RGR fNlast=fN if(nowrite.eq.0)write(8,'(5f10.5)')val,RGR,fN,SensRGR, + SensfN enddo if(nowrite.eq.0)close(8) goto 90 endif enddo write(6,'("Not a valid name; try again")') goto 101 else write(6,'("Not a valid choice. Try again")') goto 90 endif stop end * * subroutine solver(fN,RGR,fNlo0,fNhi0,nbin) implicit real*8 (a-h,o-z) real*8 Kc,Ke common alphaL,pstar,r,beta,vplant, alphaM,vM,kC,Ke,G,RGRmax,fC * Compute RGR differences for endpoints in fN 90 fNhi=fNhi0 Offhi=RGRClcalc(fNhi)-RGRNlcalc(fNhi) fNlo=fNlo0 Offlo=RGRClcalc(fNlo)-RGRNlcalc(fNlo) n=0 100 n=n+1 ! = number of halvings done to date fNmid=0.5*(fNhi+fNlo) fN=fNmid Offmid=RGRClcalc(fN)-RGRNlcalc(fN) if(Offhi*Offlo.gt.0.)then write(6,'("fNlo0, fNhi0 both give RGRCl-RGRNl of same", + " sign (",2f10.5,")",/,5x,"Root-finding likely to fail")') + Offlo,Offhi if(Offmid*Offhi.gt.0.)then write(6,'("Testing if fN range is too broad. However, there", + " is no sign change",/,5x,"in RGRCl-RGRNl even at midpoint;", + " enter now a broader range of fNlo0, fNhi0")') read(5,*)fNlo0,fNhi0 goto 90 else write(6,'("RGRCl-RGRNl changes sign twice, between fNlo0 and", + " fNhi0.",/,5x,"That is, RGRCl-RGRNl =",3f10.5," at fNlo0,", + " fNmid, fNhi0.",/,5x,"The root (correct fN) is either ", + "between fNlo0 and fNmid, or between fNmid and fNhi0.",/,5x, + "Now reset fNlo0 or fNhi0 to span the most likely interval", + " for fN.")')Offlo,Offmid,Offhi read(5,*)fNlo0,fNhi0 goto 90 endif endif if(Offhi*Offmid.lt.0.)then ! Root (true fN) lies between these two fNlo=fNmid Offlo=Offmid elseif(Offhi*Offmid.gt.0.)then fNhi=fNmid Offhi=Offmid else ! Hit root exactly (not likely!) goto 200 endif if(n.lt.nbin)goto 100 * Report final fN and other results 200 RGR=RGRClcalc(fN) write(6,300)fN,RGR 300 format("fN=",f7.5," RGR=",f7.5) PSinhib=kC/(1.+kC*fC) write(6,310)fC,PSinhib 310 format(5x,"fC=",f10.5,", PSinhib=",f8.5) if((fC.lt.0.3).or.(fC.gt.1.5))then write(6,'("fC is not in normal bounds (0.3,1.5);", + /,5x,"Is N stress to high, or are kC, Ke, RGRmax", + " unrealistic?")') endif return end * * real*8 function RGRClcalc(fN) implicit real*8 (a-h,o-z) real*8 kC,Ke common alphaL,pstar,r,beta,vplant, alphaM,vM,kC,Ke,G,RGRmax,fC fac=beta*alphaL*pstar*fN/((1.+r*(1.+G*alphaM))*RGRmax) b=1.-fac if(kC.gt.0.)then fC=(-b+dsqrt(b*b+4.*kC*fac*Ke))*0.5/kC else fC=fac*Ke/b endif RGRClcalc=beta*alphaL*pstar*fN/((1.+r*(1.+G*alphaM))*(1.+kC*fC)) return end * * real*8 function RGRNlcalc(fN) implicit real*8 (a-h,o-z) real*8 kC,Ke common alphaL,pstar,r,beta,vplant, alphaM,vM,kC,Ke,G,RGRmax,fC RGRNlcalc=r*(vplant+alphaM*vM)/((1.+r*(1.+alphaM))*fN) return end * * subroutine pstarread(pstar) implicit real*8 (a-h,o-z) character*1 opt 68 write(6,'("Options: enter letter ''A'' if you have p* in", + " desired units",/,5x,"Enter ''B'' to calc. from ", + "assim. per mass and N per leaf mass")') read(5,*,iostat=ierr)opt if(ierr.ne.0)then write(6,'("Error (probably number instead of letter);", + " try again")') goto 68 endif if(opt.eq.'A')then write(6,'("Enter pstar, in g(PSate) gN^-1 d^-1")') read(5,*)pstar elseif(opt.eq.'B')then 70 write(6,'("Enter assim. as umolCO2 g^-1 s^-1 and leaf", + " N per dry mass (e.g., 0.03 gN/gDM);",/,5x,"Also", + " enter photoperiod, in seconds")') read(5,*)Amass,fN,pper pstar=Amass*pper*30.e-6/fN 72 write(6,'("Computed pstar=",g10.5,"; should be order of", + " magnitude of about 25.",/,5x,"If correct, enter 0;", + " else, enter 1 and redo entries")')pstar read(5,*,iostat=ierr)n if(ierr.ne.0)then write(6,'("Error (probably letter instead of number);", + " try again")') goto 72 endif if(n.eq.1)goto 70 else write(6,'("Invalid option; try again (upper case A or", + " B)")') goto 68 endif return end * * subroutine vplantread(vplant) implicit real*8 (a-h,o-z) character*1 opt real*8 Km,Ja 75 write(6,'("Options: enter letter ''A'' if you have vplant", + " in desired units,",/,5x,"or ''B'' if you have it in other", + " units",/,5x,"or ''C'' to calc. from root ", + "radius, root spacing, diffusion const., etc.")') read(5,*,iostat=ierr)opt if(ierr.ne.0)then write(6,'("Error (probably number instead of letter);", + " try again")') goto 75 endif if(opt.eq.'A')then write(6,'("Enter vplant in gN gDMroot^-1 d^-1")') read(5,*)vplant elseif(opt.eq.'B')then write(6,'("Program will accept vplant in umol(N) ", + "gDMroot^-1 h^-1;",/,5x,"Please enter vplant in these ",/, + "units, and h=hours per day of active uptake")') read(5,*)vplant,h vplant=vplant*14.e-6*h 77 write(6,'("Computed vplant=",g10.5,"; should be order ", + "of magnitude of about 0.02.",/,5x,"If correct, enter", + " 0; else, enter 1 and redo entries")')vplant read(5,*,iostat=ierr)n if(ierr.ne.0)then write(6,'("Error (probably letter instead of number);", + " try again")') goto 77 endif if(n.eq.1)goto 75 elseif(opt.eq.'C')then write(6,'("We will use root radius a (m; about 0.0001), ", + "root spacing b (m; about 0.001)",/,5x, + "diffusion constant D (m^2 s^-1; about 1.e-9 for NO3-", + " in wet soil)", + /,5x,"Michaelis constant Km (mol m^-3; about 0.020 = ", + "20 uM)",/,5x,"Vmax (umolN g(DMroot)^-1 h-1, which we ", + "will convert",/,5x,"and bulk soil conc. (mol m^-3;", + " about 0.1 = 100 uM)",//,"Now enter a, b, D, Km, ", + "cb; handle Vmax separately right after this")') read(5,*)a,b,D,Km,cb 80 write(6,'("If you have Vmax in mol m^-2 s^-1, enter 1;", + /,5x," if in umol g(DMroot)^-1 h^-1, enter 2")') read(5,*,iostat=ierr)n if(ierr.ne.0)then write(6,'("Error (probably letter instead of number);", + " try again")') goto 80 endif if(n.eq.1)then write(6,'("Enter Vmax in mol m^-2 s^-1 at root", + " surface")') read(5,*)Vmax elseif(n.eq.2)then 81 write(6,'("Enter Vmax in umol g(DMroot)^-1 h^-1 AND", + " root dry mass fraction (e.g., 0.15)")') read(5,*)Vmaxugh,fDMroot Vmax=Vmaxugh*fDMroot*0.5*a/3600. * (umol g-1 h-1)*1.e-6(mol/umol)*(fDMroot gDM/gFM)*(1.e6 gDM/m3)*(a/2) * (m3 vol)/(m2 surface)*(h/s) 82 write(6,'("Computed Vmax =",g10.5," mol m-2 s-1;", + "expect about 1.e-6;",/,5x,"To accept, enter 0; to", = " redo, enter 1")')Vmax read(5,*,iostat=ierr)n if(ierr.ne.0)then write(6,'("Error (probably letter instead of number);", + " try again")') goto 82 endif if(n.eq.1)goto 81 else write(6,'("Invalid option; enter usable option")') goto 80 endif * Now calculate vplant in cylindrical diffusion geom., steady state; * first, compute steady-state conc. of NO3- at root surface, ca x=Vmax*a*dlog(b/a)/D bbig=x-cb+Km ca=(-bbig+dsqrt(bbig*bbig+4.*Km*cb))*0.5 Ja=Vmax*ca/(ca+Km) ! This is mol m-2 s-1 write(6,'(3x,"Computed ca=",f8.5," mol m-3; Ja=",6Pf9.4, + " umol m-2 s-1")')ca,Ja * Compute relative sensitivities of vplant=v to the factors Vmax, Km, * and cb, as the log. derivatives d(ln v)/d(ln Vmax), d(ln v)/d(ln Km), and * d(ln v)/d(ln cb) fac=2.*ca+bbig dKm=-Km*(ca-Km*(cb-ca)/fac)/(ca*(ca+Km)) dcb=cb*Km/(ca*fac) dVm=1.-Km*x/((ca+Km)*fac) dKm0=-Km/(cb+Km) dcb0=-dKm0 dVm0=1. write(6,'(/,5x,"Relative sensitivities of vplant to ", + "changes in Vmax, Km, or cb:",/,8x,"We express these as ", + "(relative change in vplant)/(relative change in factor)", + /,8x,"=d(ln vplant)/d(ln factor). Example: a 5% increase",/, + 8x,"in Vmax causes a 2.3% increase in vplant, which is only", + " 0.46 of 5%.",/,8x,"The sensitivity is then 0.46.",/,5x, + "For the soil and plant parameters just chosen, we have:",//, + 5x,"Sensitivity Sensitivity, IF soil With actual D:", + " expected",/,5x,"of vplant to diffus. were infinite ", + " sensitivity of RGR to",/,5x,"Vmax,Km, or cb (no ", + "depletion zones) Vmax,K, or cb (neglecting",/,47x, + "contrib. of mycorrhizae)",2x,/,"Factor",/,2x,"------",/, + 3(3x,a4,0Pf10.5,2(7x,f12.5),/),/,5x,"Note that these are", + " basically flux-control coefficients and sum to 1",/, + 8x,"over all 3 factors",/))') + "Vmax",dVm,dVm0,0.5*dVm,"Km ",dKm,dKm0,0.5*dKm, + "cb ",dcb,dcb0,0.5*dcb * We need vplant in gN g(DMroot)^-1 d^-1; convert if(n.eq.1)then ! We did not previously enter fDMroot write(6,'("Enter root dry mass fraction (e.g., 0.15)")') read(5,*)fDMroot endif write(6,'("Enter length of uptake period per day, in h ", + "(say, 14)")') read(5,*)pper vplant=Ja*14.*(2./a)*(1.e-6/fDMroot)*pper*3600. 84 write(6,'("Computed vplant=",g10.5,"; should be order ", + "of magnitude of about 0.02.",/,5x,"If correct, enter", + " 0; else, enter 1 and redo entries")')vplant read(5,*,iostat=ierr)n if(ierr.ne.0)then write(6,'("Error (probably letter instead of number);", + " try again")') goto 84 endif if(n.eq.1)goto 75 else write(6,'("Invalid option; try again (upper case A or", + " B)")') goto 75 endif ! End of block-if to handle the 2 options for vplant return end * * subroutine vMread(vplant,vM) implicit real*8 (a-h,o-z) character*1 opt 68 write(6,'("Options: enter letter ''A'' if you have vM in", + " desired units",/,5x,"Enter ''B'' to calc. as vplant", + " multiplied by a factor")') read(5,*,iostat=ierr)opt if(ierr.ne.0)then write(6,'("Error (probably number instead of letter);", + " try again")') goto 68 endif if(opt.eq.'A')then write(6,'("Enter vM, in gN g(DMroot)^-1 d^-1")') read(5,*)vM elseif(opt.eq.'B')then 70 write(6,'("Enter multiplier, vM/vplant")') read(5,*)vmult vM=vplant*vmult else write(6,'("Invalid option; try again (upper case A or", + " B)")') goto 68 endif return end * * * Perhaps later: make option 'C' for vM: determine vM from diffusion * kinetics and Vmax, Km of MR fungi. Use the subroutine vread(vM), * which will be vplantread generalized to read either vplant or vM. * This generalization includes changing the prompts to read, for * example, per mass of hyphae or per mass of root, as appropriate * ** To do later: Use pointers (f90 style) instead of equivalence stmt. in ** subroutine paramread ** Later - it's complicated, and not all people will have f90 compilers