!################################################################# ! ! routine to perform main processing and ! handles time ordering of model processes ! !### macro's ##################################################### ! #define TRACEBACK write (gol,'("in ",a," (",a,", line",i5,")")') rname, __FILE__, __LINE__; call goErr #define IF_NOTOK_RETURN(action) if (status/=0) then; TRACEBACK; action; return; end if #define IF_ERROR_RETURN(action) if (status> 0) then; TRACEBACK; action; return; end if ! #include "tm5.inc" ! !################################################################# module modelIntegration use GO, only : gol, goPr, goErr, goBug, goLabel implicit none ! --- in/out ----------------------------- private public :: Proces_Init, Proces_Done, Proces_Update public :: Proces_Region ! local .. character(len=1) :: output_after_step ! --- const ---------------------------------- ! module name character(len=*), parameter :: mname = 'ModelIntegration' ! --- var ---------------------------------- ! timer handles: integer :: itim_advectx, itim_advecty, itim_advectz integer :: itim_vertical, itim_chemistry, itim_source_sink contains ! ===================================================== subroutine Proces_Init( status ) use GO, only : GO_Timer_Def use GO, only : TrcFile, Init, Done, ReadRc use global_data, only : rcfile use Binas , only : p_global use dims , only : nregions, nregions_all, iglbsfc use global_data, only : rcfile use Meteo , only : sp_dat, phlb_dat, m_dat use Meteo , only : sp1_dat, sp2_dat, tsp_dat, spm_dat use Meteo , only : mfu_dat, mfv_dat use Meteo , only : mfw_dat use Meteo , only : temper_dat, humid_dat use Meteo , only : entu_dat, entd_dat, detu_dat, detd_dat use Meteo , only : lwc_dat, iwc_dat, cc_dat, cco_dat, ccu_dat use Meteo , only : gph_dat, omega_dat use Meteo , only : oro_dat, lsmask_dat use Meteo , only : Set use Meteo , only : Meteo_Alloc #ifndef without_advection use Advect , only : Advect_Init #endif #ifndef without_convection use Convection , only : Convection_Init #ifndef without_diffusion use Diffusion , only : Diffusion_Init #endif #endif #ifndef without_chemistry use Chemistry , only : Chemistry_Init #endif use Sources_Sinks , only : Sources_Sinks_Init #ifndef without_dry_deposition use Dry_Deposition, only : Dry_Deposition_Init #endif #ifndef without_wet_deposition use Wet_Deposition, only : Wet_Deposition_Init #endif ! --- in/out -------------------------------- integer, intent(out) :: status ! --- const ------------------------------- character(len=*), parameter :: rname = mname//'/Proces_Init' ! --- local --------------------------------- type(TrcFile) :: rcF integer :: n ! --- begin --------------------------------- call goLabel(rname) ! ! select meteo to be used ! ! loop over all (zoom) regions: do n = 1, nregions ! current (surface) pressure and mass: call Set( sp_dat(n), status, used=.true. ) call Set( m_dat(n), status, used=.true. ) call Set( phlb_dat(n), status, used=.true. ) ! surface pressures, tendency, and average (mid) call Set( sp1_dat(n), status, used=.true. ) call Set( sp2_dat(n), status, used=.true. ) call Set( tsp_dat(n), status, used=.true. ) call Set( spm_dat(n), status, used=.true. ) end do ! ! check fpp macro's ! ! macro 'with_zoom' should be defined if nregions>1 if ( nregions > 1 ) then #ifndef with_zoom write (gol,'("macro with_zoom not defined while nregons = ",i2)') nregions; call goErr TRACEBACK; status=1; return #endif end if ! chance to have faster code ... if ( nregions == 1 ) then #ifdef with_zoom write (gol,'("WARNING - macro with_zoom defined while nregions=1; slower run")'); call goPr #endif end if ! ! init processes ! ! NOTE: these inits should select the appropriate meteo #ifndef without_advection ! init advection ? call Advect_Init( status ) IF_NOTOK_RETURN(status=1) #endif #ifndef without_convection ! init convec: call Convection_Init( status ) IF_NOTOK_RETURN(status=1) #ifndef without_diffusion ! init diffusion: call Diffusion_Init( status ) IF_NOTOK_RETURN(status=1) #endif #endif #ifndef without_chemistry ! init chemistry: call Chemistry_Init( status ) IF_NOTOK_RETURN(status=1) #endif ! init sources and sinks: call Sources_Sinks_Init( status ) IF_NOTOK_RETURN(status=1) #ifndef without_dry_deposition ! setup dry depos: call Dry_Deposition_Init( status ) IF_NOTOK_RETURN(status=1) #endif #ifndef without_wet_deposition ! setup wet depos: call Wet_Deposition_Init( status ) IF_NOTOK_RETURN(status=1) #endif ! ! allocate meteo ! ! some meteo implies other ... do n = 1, nregions_all ! convec requires gph and omega: if ( entu_dat(n)%used ) then call Set( gph_dat(n), status, used=.true. ) call Set( omega_dat(n), status, used=.true. ) end if ! omega (Pa/s) is computed form mfw (kg/s): if ( omega_dat(n)%used ) then call Set( mfw_dat(n), status, used=.true. ) end if ! gph is computed from oro/sp/temper/humid if ( gph_dat(n)%used ) then call Set( oro_dat(n), status, used=.true. ) call Set( sp_dat(n), status, used=.true. ) call Set( temper_dat(n), status, used=.true. ) call Set( humid_dat(n), status, used=.true. ) end if ! sp is interpolated in time between sp1 and sp2: if ( sp_dat(n)%used ) then call Set( sp1_dat(n), status, used=.true. ) call Set( sp2_dat(n), status, used=.true. ) end if ! cloud covers and overhead/underfeet cloud covers if ( cc_dat(n)%used ) then call Set( cco_dat(n), status, used=.true. ) call Set( ccu_dat(n), status, used=.true. ) end if end do !>>> Please use obscure meteo fields only if you really need them ... ! If too many fields are used by default, even a simple meteo production run ! of an extra surface field will request fields that are not used anyway ... ! write (gol,'("WARNING - ")'); call goPr write (gol,'("WARNING - global orography not in use by default anymore;")'); call goPr write (gol,'("WARNING - see the comment in ",a)') rname; call goPr write (gol,'("WARNING - ")'); call goPr ! ! always store also the orography on 1x1: ! call Set( oro_dat(iglbsfc), status, used=.true. ) ! !<<< ! allocate used meteo call Meteo_Alloc( status ) IF_NOTOK_RETURN(status=1) ! ! default values ! ! Fill default values for surface pressure; ! might be required for mass fields used in init routines do n = 1, nregions_all if ( sp_dat(n)%used ) sp_dat(n)%data = p_global if ( spm_dat(n)%used ) spm_dat(n)%data = p_global end do ! ! output sampling ! ! read how to sample the output: ! x,y,z,c,v,e : always after this process ! a : after all steps (testing, not recommended) ! d : the 'old' way (default) call Init( rcF, rcfile, status ) IF_NOTOK_RETURN(status=1) call ReadRc( rcF, 'output.after.step', output_after_step, status, default = 'd' ) IF_ERROR_RETURN(status=1) if (status == 0) then write(gol, *) '****************************************' ; call goPr write(gol, *) 'Output will be collected after step : ',output_after_step ; call goPr write(gol, *) '****************************************' ; call goPr elseif ( status < 0) then write(gol, *) '**********************************************' ; call goPr write(gol, *) 'Output will be collected in the OLD way' ; call goPr write(gol, *) 'consider to include e.g. output.after.step : v' ; call goPr write(gol, *) '**********************************************' ; call goPr endif ! close rcfile: call Done( rcF, status ) IF_NOTOK_RETURN(status=1) ! ! timing ! ! define ... call GO_Timer_Def( itim_advectx , 'advectx' , status ) IF_NOTOK_RETURN(status=1) call GO_Timer_Def( itim_advecty , 'advecty' , status ) IF_NOTOK_RETURN(status=1) call GO_Timer_Def( itim_advectz , 'advectz' , status ) IF_NOTOK_RETURN(status=1) call GO_Timer_Def( itim_vertical , 'vertical' , status ) IF_NOTOK_RETURN(status=1) call GO_Timer_Def( itim_chemistry , 'chemistry' , status ) IF_NOTOK_RETURN(status=1) call GO_Timer_Def( itim_source_sink, 'source_sink' , status ) IF_NOTOK_RETURN(status=1) ! ! done ! ! end call goLabel() ! ok status = 0 end subroutine Proces_Init ! *** subroutine Proces_Done( status ) #ifndef without_chemistry use Chemistry , only : Chemistry_Done #endif #ifndef without_advection use Advect , only : Advect_Done #endif #ifndef without_convection use Convection , only : Convection_Done #ifndef without_diffusion use Diffusion , only : Diffusion_Done #endif #endif use Sources_Sinks , only : Sources_Sinks_Done #ifndef without_dry_deposition use Dry_Deposition, only : Dry_Deposition_Done #endif #ifndef without_wet_deposition use Wet_Deposition, only : Wet_Deposition_Done #endif ! --- in/out -------------------------------- integer, intent(out) :: status ! --- const ------------------------------- character(len=*), parameter :: rname = mname//'/Proces_Done' ! --- begin --------------------------------- #ifndef without_advection call Advect_Done( status ) IF_NOTOK_RETURN(status=1) #endif #ifndef without_convection ! done with convec: call Convection_Done( status ) IF_NOTOK_RETURN(status=1) #ifndef without_diffusion ! done with diffusion: call Diffusion_Done( status ) IF_NOTOK_RETURN(status=1) #endif #endif #ifndef without_chemistry ! done with chemistry: call Chemistry_Done( status ) IF_NOTOK_RETURN(status=1) #endif ! done with sources and sinks: call Sources_Sinks_Done( status ) IF_NOTOK_RETURN(status=1) #ifndef without_dry_deposition ! done with dry depos call Dry_Deposition_Done( status ) IF_NOTOK_RETURN(status=1) #endif #ifndef without_wet_deposition ! done with wet depos call Wet_Deposition_Done( status ) IF_NOTOK_RETURN(status=1) #endif ! ok status = 0 end subroutine Proces_Done ! *** ! ! update process fields, for example because meteo has been changed ! subroutine Proces_Update( status ) use dims , only : nregions #ifndef without_diffusion use Diffusion , only : Calc_Kzz, diffusion_write, read_diffusion, write_diffusion use toolbox , only : escape_tm #endif use sources_sinks , only : trace_after_read #ifndef without_dry_deposition use dry_deposition, only : dd_surface_fields #endif #ifndef without_wet_deposition use wet_deposition, only : calculate_lsp_scav #endif ! --- in/out -------------------------------- integer, intent(out) :: status ! --- const ------------------------------- character(len=*), parameter :: rname = mname//'/Proces_Update' ! --- local -------------------------------- integer :: region ! --- begin --------------------------------- #ifndef without_diffusion ! recalculate diffusion coeffiecents if necessary ! If diffusion_write, then read existing DKG file or ! write it if it does not already exist. if (diffusion_write) then call read_diffusion( status ) if ( status > 0 ) then ! error reading call escape_tm('Tracer/calc_kzz: error reading diffusion') else if ( status < 0 ) then ! dkg file could not be found --> calculate call Calc_Kzz ! write to file call write_diffusion( status ) if ( status /= 0 ) call escape_tm('Tracer: error writing diffusion') end if else call Calc_Kzz endif #endif #ifndef without_dry_deposition ! calculate the dry_deposition fields from stuff in fsurf file (will coarsen...) call dd_surface_fields( status ) IF_NOTOK_RETURN(status=1) #endif #ifndef without_wet_deposition ! pre-calculate the lsp scavenging rloss1, rloss2, rloss3 do region = 1, nregions call calculate_lsp_scav(region) end do #endif ! recalculate chemistry fields if necessary call trace_after_read( status ) IF_NOTOK_RETURN(status=1) ! ok status = 0 end subroutine Proces_Update ! *** ! ! for a given region, performs tref(region)/tref(parent(region)) ! advection steps for each direction x,y,z ! including all region's children, granchildren, etc. ! written by patrick berkvens and mike botchev, march-june 1999 !WP! Also chemistry and sources_sinks are called during advection, !WP! wouter peters, july 2000 ! !example (cmk) ! XYZ ECCEZYX ! xyz eccezyx cezyx xyz ec ! xyzeccezyx ecxyzzyxec cezyxxyzec zyxeccexyz ! in this case the operations e and c should be executed in the ! interface region ! but not in the core of the zoom region (otherwise double counting) ! This results in the most simple algorithm, because you should leave ! the DO_STEPS routine AFTER ! (1) finishing all the steps... ! (2) XYZ ! MK 2000... ! recursive subroutine Proces_Region( region, tr, status ) use GO , only : TDate, NewDate use dims , only : tref, revert, ndyn, itaur use dims , only : parent, children use dims , only : statusr => status use dims , only : n_operators use datetime , only : tau2date use user_output , only : user_output_step use zoom_tools , only : update_parent use advect_tools, only : m2phlb, m2phlb1 use tracer_data , only : Par_Check_Domain ! --- in/out ------------------------------- integer, intent(in) :: region type(TDate) :: tr(2) integer, intent(out) :: status ! --- const ------------------------------- character(len=*), parameter :: rname = mname//'/proces_region' ! --- local --------------------------------- integer :: child, i, ichild, tref_, n_children, dtime integer :: time6(6) type(TDate) :: tri(2) ! --- begin ---------------------------------- ! determine refinement factor with respect to the parent tref_ = tref(region)/tref(parent(region)) !main time step for all processes. note that all timesteps are default/2 dtime = revert*ndyn/(2*tref(region)) !e.g. xyzvec (firts call) and cevzyx (second call) do i = 1, tref_ ! time step: itaur(region)+[0,dtime] call tau2date( itaur(region), time6 ) tri(1) = NewDate( time6=time6 ) call tau2date( itaur(region)+dtime, time6 ) tri(2) = NewDate( time6=time6 ) call do_steps( region, tri, status ) IF_NOTOK_RETURN(status=1) ! CALL advect_region for all the children (IF there are any) ! ! children(region,0 ) : number of children ! children(region,1:children(region,0)) : list with region numbers ! ! example for glb6x4/eur3x2/eur1x1 : ! children(1,:)= 1, 2 0 0 ! children(2,:)= 1, 3 0 0 ! children(3,:)= 0, 0 0 0 ! ichild = 0 do while ( ichild < children(region,0) ) ichild = ichild + 1 child = children(region,ichild) call proces_region( child, tri, status ) IF_NOTOK_RETURN(status=1) end do ! do the remaining steps if necessary...! if ( mod(statusr(region),n_operators) /= 0 ) then call do_steps( region, tri, status ) IF_NOTOK_RETURN(status=1) end if ! update region time..... itaur(region) = itaur(region) + dtime !>>> moved back here (MK) !>>> also moved to de_steps for specific sampling..... !! Accumulate output, write info at stations !! Note: the edges in the interface may not have been updated: !! this is the task of the parent. if (output_after_step == 'd') then !old default call Par_Check_domain(region,'n','tracer') call user_output_step( region, status ) IF_NOTOK_RETURN(status=1) endif !<<< end do ! update parent using child data ? ! compute pressure grid from (changed) air mass: if ( region /= 1 ) then call update_parent( region ) call m2phlb( region ) else call m2phlb1( 1 ) end if !if(okdebug) call Par_Check_mass(region,'after proces_region') ! ok status = 0 end subroutine proces_region ! *** !------------------------------------------------------------- ! performs next three (x-, y- and z-) advection steps ! written by patrick berkvens and mike botchev, march-june 1999 !wp! added the call to chemistry. subroutine now does four !wp! steps (c-x-y-z) wp, july 2000 !mk! changed the order for proper processing of the interfaces. !mk! feb-2001 !------------------------------------------------------------- subroutine do_steps( region, tr, status ) use GO , only : TDate use GO , only : GO_Timer_Start, GO_Timer_End use dims , only : statusr => status use dims , only : okdebug use dims , only : splitorderzoom, n_operators #ifndef without_chemistry use chemistry , only : chemie use sources_sinks , only : source2 #endif use sources_sinks , only : source1 #ifdef with_budgets use budget_global , only : budget_transportg #endif use toolbox , only : escape_tm #ifdef MPI use mpi_const , only : previous_par #endif use ParTools , only : myid, root use ParTools , only : ntracetloc, ntracet_ar use tracer_data , only : Par_Check_Domain use global_data , only : mass_dat use chem_param , only : names, ntracet #ifndef without_convection use convection , only : convec #endif #ifndef without_advection use advectx , only : advectxzoom, put_xedges use advecty , only : advectyzoom, put_yedges use advectz , only : advectzzoom #endif #ifndef without_wet_deposition use wet_deposition, only : lspscav #endif #ifdef with_tendencies use tracer_data, only : PLC_Set, plc_ntr, plc_ipr_tchem #endif use user_output , only : user_output_step ! --- in/out ----------------------------------- integer, intent(in) :: region type(TDate) :: tr(2) integer, intent(out) :: status ! --- const ------------------------------- character(len=*), parameter :: rname = mname//'/do_steps' ! --- local ----------------------------------- integer :: child, i123, ichild, tref_child, reg, rgi character :: tobedone character(len=1) :: next_step, prev_step !cmk integer :: start_x, stop_x, start_y, stop_y, n !CMKTEMP #ifdef with_tendencies integer :: itr #endif ! --- begin ----------------------------------- if ( okdebug .and. myid==root ) print *,'do_steps: region ',region !example (cmk) ! XYZ EC CEX ! xyzeccezyx cexyzzyxec ! in this case the operations e and c should be executed in ! the interface region but not in the core of the zoom region ! (otherwise double counting) ! This results in the most simple algorithm, because you should leave ! the DO_STEPS routine AFTER a Z OR after finishing all steps... ! note that the parent is responsible for the interface with the children ! this means that the interface may not be 'update' in case of (e.g.) ! xyz a ! xyzaazyx ! in this case, the edges of the child are not updated with the proces a. ! this has consequences for the mmix and save output.... ! THE INTERFACE IS PART OF THE PARENT, NOT THE CHILD..... do i123=1,n_operators !WP! changed 3 to n_operators next_step = splitorderzoom(region,statusr(region)+1) if ( statusr(region) /= 0 ) then prev_step = splitorderzoom(region,statusr(region)) else prev_step = ' ' end if if ( okdebug .and. myid == root ) then ! it's only to make work of the code visible -- step to be done ! will be printed with capital letter (X,Y or Z) do reg=1,region tobedone = ' ' if ( reg == region ) tobedone = upper(next_step) print *, 'do_steps: ',reg,': ', & splitorderzoom(reg,1:statusr(reg)),tobedone end do end if select case(next_step) case ( 'v' ) call GO_Timer_Start( itim_vertical, status ) IF_NOTOK_RETURN(status=1) #ifndef without_convection !WP! data must be on tracers for convection call Par_Check_domain(region,'n','tracer') !if ( okdebug ) call Par_Check_mass(region,'bef_convec') #ifdef with_budgets call budget_transportg(region,0,'conv',prev_step) #endif !if ( okdebug ) print *, 'do_steps: Call convection' call convec(region, status ) IF_NOTOK_RETURN(status=1) #ifdef with_budgets call budget_transportg(region,1,'conv',prev_step) #endif !if ( okdebug ) call Par_Check_mass(region,'aft_convec') #endif call GO_Timer_End( itim_vertical, status ) IF_NOTOK_RETURN(status=1) case ( 'c' ) call GO_Timer_Start( itim_chemistry, status ) IF_NOTOK_RETURN(status=1) #ifndef without_chemistry !WP! data must be on levels for chemistry call Par_Check_domain(region,'n','levels') !if ( okdebug ) call Par_Check_mass(region,'bef_chem') #ifdef with_tendencies ! store concentrations before chemistry in plc process '0' do itr = 1, plc_ntr call PLC_Set( 'fill-tracer', itr, 0, status ) IF_NOTOK_RETURN(status=1) end do #endif if ( okdebug ) print*, 'do_steps: Call chemie' call Chemie( region, tr, status ) IF_NOTOK_RETURN(status=1) !if ( okdebug ) call Par_Check_mass(region,'aft_chem') #ifdef with_tendencies ! add change due to chemistry to 'tchem' do itr = 1, plc_ntr call PLC_Set( 'add-change', itr, plc_ipr_tchem, status, ipr2=0, t=tr(2) ) IF_NOTOK_RETURN(status=1) end do #endif #endif call GO_Timer_End( itim_chemistry, status ) IF_NOTOK_RETURN(status=1) case( 's' ) call GO_Timer_Start( itim_source_sink, status ) IF_NOTOK_RETURN(status=1) !WP! data must be on tracers for convection call Par_Check_domain(region,'n','tracer') !if ( okdebug ) call Par_Check_mass(region,'bef_source') if ( okdebug ) print*, 'do_steps: Call source' call source1( region, tr, status ) IF_NOTOK_RETURN(status=1) #ifndef without_chemistry call source2( region, tr, status ) IF_NOTOK_RETURN(status=1) #endif #ifndef without_wet_deposition ! wet removal by lsp; note that dry deposition is applied in chemistry call lspscav( region ) #endif !if ( okdebug ) call Par_Check_mass(region,'aft_source') call GO_Timer_End( itim_source_sink, status ) IF_NOTOK_RETURN(status=1) case( 'x' ) call GO_Timer_Start( itim_advectx, status ) IF_NOTOK_RETURN(status=1) #ifndef without_advection !WP! data must be on tracers for X call Par_Check_domain(region,'c','tracer') !if ( okdebug ) call Par_Check_mass(region,'befX_tracers') #ifdef with_budgets call budget_transportg(region,0,'advx',prev_step) #endif call advectxzoom(region) #ifdef with_budgets call budget_transportg(region,1,'advx',prev_step) #endif !if ( okdebug ) call Par_Check_mass(region,'aftX_tracers') #endif call GO_Timer_End( itim_advectx, status ) IF_NOTOK_RETURN(status=1) case( 'y' ) call GO_Timer_Start( itim_advecty, status ) IF_NOTOK_RETURN(status=1) #ifndef without_advection !WP! data must be on tracers for Y call Par_Check_domain(region,'c','tracer') !if ( okdebug ) call Par_Check_mass(region,'befY_tracers') #ifdef with_budgets call budget_transportg(region,0,'advy',prev_step) #endif call advectyzoom(region) #ifdef with_budgets call budget_transportg(region,1,'advy',prev_step) #endif !if ( okdebug ) call Par_Check_mass(region,'aftY_tracers') #endif call GO_Timer_End( itim_advecty, status ) IF_NOTOK_RETURN(status=1) case( 'z' ) call GO_Timer_Start( itim_advectz, status ) IF_NOTOK_RETURN(status=1) #ifndef without_advection !WP! data must be on tracers for Z call Par_Check_domain(region,'n','tracer') !if ( okdebug ) call Par_Check_mass(region,'befZ_tracers') #ifdef with_budgets call budget_transportg(region,0,'advz',prev_step) #endif call advectzzoom(region) #ifdef with_budgets call budget_transportg(region,1,'advz',prev_step) #endif !if ( okdebug ) call Par_Check_mass(region,'aftZ_tracers') #endif call GO_Timer_End( itim_advectz, status ) IF_NOTOK_RETURN(status=1) case default print *,'do_steps: strange value in splitorderzoom: ', next_step print *,'do_steps: (must be c, e, v, x, y or z)' call escape_tm('do_steps: Error') end select statusr(region) = statusr(region)+1 ! e.g. xyzvec...........cevzyx ! exit after xyz or vec of cevzyx.... ! if the exit is after vec...the interface of child should be updated... ! if ( mod(statusr(region),n_operators) == 0 ) then ! after xyz....vec...update child interface CMK removed feb 2003... ! BUT replaced in march 2003 by MK: budgets were corrupted when ! the edges are NOT updated!!! still to find out why??? #ifndef without_advection if ( next_step /= 'x' ) then call put_xedges(region) call put_yedges(region) end if ! if not .... zyxcev #endif exit ! e.g. vceecvzyx end if ! flexible sample scheme ('a' means after all steps), d = old 'default', ! output_after_step = v,c,e,x,y,z,d,a if ((output_after_step == next_step) .or. (output_after_step == 'a')) then call Par_Check_domain(region,'n','tracer') call user_output_step( region, status ) IF_NOTOK_RETURN(status=1) endif !exit after 'yz' sequence to proces xyz children if ( next_step == 'z' ) then prev_step = splitorderzoom(region,statusr(region)-1) if ( prev_step == 'y' ) exit end if end do ! ok status = 0 end subroutine do_steps ! *** character function upper(xyz) implicit none character(1),intent(in) :: xyz if (xyz=='x') then upper = 'X' else if (xyz=='y') then upper = 'Y' else if (xyz=='z') then upper = 'Z' else if (xyz=='c') then upper = 'C' else if (xyz=='s') then upper = 'S' else if (xyz=='v') then upper = 'V' else upper = '_' end if end function upper end module modelIntegration