#include "forttype.h" module mod_spec_level1 use mod_parameters use mod_template use mod_process use mod_mem_alloc use mod_lunits implicit none contains !***************************************************************** subroutine add_template_spec_level1 implicit none type (template), pointer :: add character, dimension(:), allocatable :: sizestring __INTEGER :: ierr, ierr1, ierr2, ierr3, ierr4, ierr5 allocate( add, stat=ierr ) if (ierr /= 0) then write(lu_log,*) " add_template_spec_level1 : allocation failed" call peg_error('bad allocate stat') end if #ifdef MEMORY_TRACKING call mem_track('add_template_spec_level1','add',sizeoftemplate) #endif add%name = 'spec_level1' add%num_inputs = 9 allocate( add%input(add%num_inputs), stat=ierr ) if (ierr /= 0) then write(lu_log,*) " add_template_spec_level1 : allocation failed" call peg_error('bad allocate stat') end if #ifdef MEMORY_TRACKING call mem_track('add_template_spec_level1','add%input',add%num_inputs*sizeoftemplateinput) #endif allocate( add%input_type(add%num_inputs), stat=ierr ) if (ierr /= 0) then write(lu_log,*) " add_template_spec_level1 : allocation failed" call peg_error('bad allocate stat') end if #ifdef MEMORY_TRACKING call mem_track('add_template_spec_level1','add%input_type',add%num_inputs*sizeoftemplateinput_type) #endif add%input(1) = 'spec_level1_inp' add%input_type(1) = 1 add%input(2) = 'int1' add%input_type(2) = 1 add%input(3) = 'x' add%input_type(3) = 1 add%input(4) = 'int1' add%input_type(4) = 5 add%input(5) = 'x' add%input_type(5) = 5 add%input(6) = 'hole' add%input_type(6) = 5 add%input(7) = 'iregion_inp' add%input_type(7) = 5 add%input(8) = 'intrp' add%input_type(8) = 6 add%input(9) = 'level1_list' add%input_type(9) = 1 add%num_outputs = 2 allocate( add%output(add%num_outputs), stat=ierr ) if (ierr /= 0) then write(lu_log,*) " add_template_spec_level1 : allocation failed" call peg_error('bad allocate stat') end if #ifdef MEMORY_TRACKING call mem_track('add_template_spec_level1','add%output',add%num_outputs*sizeoftemplateoutput) #endif allocate( add%output_type(add%num_outputs), stat=ierr ) if (ierr /= 0) then write(lu_log,*) " add_template_spec_level1 : allocation failed" call peg_error('bad allocate stat') end if #ifdef MEMORY_TRACKING call mem_track('add_template_spec_level1','add%output_type',add%num_outputs*sizeoftemplateoutput_type) #endif add%output(1) = 'level1' add%output_type(1) = 6 add%output(2) = 'orphan' add%output_type(2) = 1 call add_template( add ) call add_template( add ) deallocate( add%input, stat=ierr1 ) deallocate( add%input_type, stat=ierr2 ) deallocate( add%output, stat=ierr3 ) deallocate( add%output_type, stat=ierr4 ) deallocate( add, stat=ierr5 ) if ( ierr1.ne.0 .or. ierr2.ne.0 .or. ierr3.ne.0 .or. & ierr4.ne.0 .or. ierr5.ne.0 ) then call peg_error('bad deallocate stat : add_template_spec_level1') endif #ifdef MEMORY_TRACKING call mem_track('add_template_spec_level1','add%input', -1) call mem_track('add_template_spec_level1','add%input_type', -1) call mem_track('add_template_spec_level1','add%output', -1) call mem_track('add_template_spec_level1','add%output_type', -1) call mem_track('add_template_spec_level1','add', -1) #endif return end subroutine add_template_spec_level1 !***************************************************************** subroutine add_process_spec_level1 use mod_parameters use mod_mesh use mod_xvals use mod_spec_level1_inp use mod_mem_alloc implicit none type (process), pointer :: add __INTEGER :: m1, m2, n __INTEGER :: ierr, ierr1, ierr2, ierr3 character (len=nchar2) :: file_spec_level1_inp logical :: exist_spec_level1_inp, write_spec_level1_inp __REAL :: qtol_tmp, qcutoff_tmp character, dimension(:), allocatable :: sizestring do m1 = 1,nmesh if( phantom(m1) ) cycle allocate( add, stat=ierr ) if (ierr /= 0) then write(lu_log,*) " add_process_spec_level1 : allocation failed" call peg_error('bad allocate stat') end if #ifdef MEMORY_TRACKING call mem_track('add_process_spec_level1','add', sizeofprocess) #endif add%name = 'spec_level1' add%primary_1 = name(m1) add%primary_2 = blank_string n = 0 do m2 = 1,nmesh if( phantom(m2) ) cycle if( check_intersection( m1, m2 ) .and. m1.ne.m2 ) n = n + 1 enddo add%num_secondary_1 = n allocate( add%secondary_1(add%num_secondary_1), stat=ierr ) if (ierr /= 0) then write(lu_log,*) " add_process_spec_level1 : allocation failed" call peg_error('bad allocate stat') end if #ifdef MEMORY_TRACKING call mem_track('add_process_spec_level1','add%secondary_1',add%num_secondary_1*sizeofprocesssecondary_1) #endif n = 0 do m2 = 1,nmesh if( phantom(m2) ) cycle if( check_intersection( m1, m2 ) .and. m1.ne.m2 ) then n = n + 1 add%secondary_1(n) = name(m2) endif enddo call check_level1_list( add ) add%num_secondary_2 = 0 allocate( add%secondary_2(add%num_secondary_2), stat=ierr ) if (ierr /= 0) then write(lu_log,*) " add_process_spec_level1 : allocation failed" call peg_error('bad allocate stat') end if #ifdef MEMORY_TRACKING call mem_track('add_process_spec_level1','add%secondary_2',add%num_secondary_2*sizeofprocesssecondary_2) #endif call add_process( add ) deallocate( add%secondary_1, stat=ierr1) deallocate( add%secondary_2, stat=ierr2) deallocate( add, stat=ierr3) if ( ierr1.ne.0 .or. ierr2.ne.0 .or. ierr3.ne.0 ) then call peg_error('bad deallocate stat : add_template_spec_level1') endif #ifdef MEMORY_TRACKING call mem_track('add_process_spec_level1','add%secondary_1',-1) call mem_track('add_process_spec_level1','add%secondary_2',-1) call mem_track('add_process_spec_level1','add',-1) #endif enddo do m1 = 1, nmesh if( phantom(m1) ) cycle write_spec_level1_inp = .false. call set_file_name_spec_level1_inp( name(m1), file_spec_level1_inp ) inquire( file=file_spec_level1_inp, exist=exist_spec_level1_inp ) if( exist_spec_level1_inp ) then call load_spec_level1_inp( name(m1), qtol_tmp, qcutoff_tmp ) if( qtol(m1) .ne. qtol_tmp .or. & qcutoff(m1) .ne. qcutoff_tmp ) write_spec_level1_inp = .true. else write_spec_level1_inp = .true. endif if( write_spec_level1_inp ) then call store_spec_level1_inp( name(m1), qtol(m1), qcutoff(m1) ) endif enddo return end subroutine add_process_spec_level1 !***************************************************************** subroutine check_level1_list( add ) use mod_parameters use mod_mesh use mod_lunits use mod_mem_alloc implicit none type (process), pointer :: add __INTEGER :: date, time __INTEGER :: m, n, nh, num_level1_tmp logical :: exist_level1_list, write_level1_list, match character (len=nchar2) :: file_level1_list character (len=nchar) :: mesh_tmp m = mesh_no(add%primary_1) write_level1_list = .false. call set_file_name_level1_list( name(m), file_level1_list ) inquire( file=file_level1_list, exist=exist_level1_list ) if( exist_level1_list ) then open( lu_level1_list, file=file_level1_list, status='old', form='unformatted' ) read(lu_level1_list) read(lu_level1_list)num_level1_tmp if( add%num_secondary_1 .eq. num_level1_tmp ) then do n = 1,num_level1_tmp read(lu_level1_list)mesh_tmp match = .false. do nh = 1,add%num_secondary_1 if( add%secondary_1(nh) .eq. mesh_tmp) then match = .true. exit endif enddo if( .not. match ) then write_level1_list = .true. exit endif enddo else write_level1_list = .true. endif close( lu_level1_list ) else write_level1_list = .true. endif if( write_level1_list ) then open( lu_level1_list, file=file_level1_list, status='unknown', form='unformatted' ) call get_time( date, time ) write(lu_level1_list) date, time write(lu_level1_list) add%num_secondary_1 do n = 1,add%num_secondary_1 write(lu_level1_list) add%secondary_1(n) enddo write(lu_level1_list) -date, -time close(lu_level1_list) endif return end subroutine check_level1_list !***************************************************************** subroutine set_file_name_level1_list( name, file_name ) use mod_dir use mod_mem_alloc implicit none character(len=*) :: name character(len=*) :: file_name file_name = trim(level1_dir)//'/'//trim(name)//'.level1_list' return end subroutine set_file_name_level1_list !***************************************************************** subroutine execute_spec_level1( proc ) use mod_ptr_x use mod_x use mod_ptr_intrp use mod_mesh use mod_intrp use mod_xvals use mod_lunits use mod_log_file use mod_ptr_int1 use mod_ptr_hole use mod_hole use mod_int1 use mod_level1 use mod_orphan use mod_ptr_project use mod_project use mod_shift use mod_mem_alloc use mod_region use mod_volume implicit none type (process), pointer :: proc __INTEGER, dimension(2) :: jg, kg, lg __INTEGER, dimension(:), allocatable :: iint1, nfr, msh, ist __INTEGER, dimension(:,:,:), allocatable :: iblank, iblank_tmp __INTEGER :: m1, m2, ierr, m, n, nl, nr, nv __INTEGER :: jb, kb, lb, ji, ki, li __INTEGER :: j,k,l __INTEGER :: num_iint1, num_iint1f1, num_iint1f2, num_iint1f1b, num_iint1f2b __INTEGER :: num_isten, num_iorphan, num_project, num_ihole, num_ilevel1 __INTEGER :: jbeg, jend, kbeg, kend, lbeg, lend __REAL, dimension(:), allocatable :: ql1, ql2, cd, dxt, dyt, dzt __REAL :: qual1_tmp, qual2_tmp, cell_diff_tmp __REAL :: dxint_tmp, dyint_tmp, dzint_tmp __REAL :: xp, yp, zp, cdr __REAL :: time, date type (level1) :: l1data type (orphan) :: orph_data logical :: intest, stntest m1 = mesh_no(proc%primary_1) allocate( iblank_tmp(mjmax(m1),mkmax(m1),mlmax(m1)), stat=ierr) if( ierr /= 0 ) then write(lu_log,*) " iblank: allocation failed" call peg_error('bad allocate stat') endif #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','iblank_tmp',mjmax(m1)*mkmax(m1)*mlmax(m1)*sizeofint) #endif iblank_tmp=1 ! ! Read in lists of first and second layer fringe points ! call load_int1( proc%primary_1, num_iint1f1, iint1f1, & num_iint1f2, iint1f2 ) num_iint1 = num_iint1f1 + num_iint1f2 ! ! Allocate arrays to hold fringe point indicies ! allocate( iint1(num_iint1), nfr(num_iint1), stat=ierr) if( ierr /= 0 ) then write(lu_log,*) " iint1: allocation failed" call peg_error('bad allocate stat') endif #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','iint1',num_iint1*sizeofint) call mem_track('execute_spec_level1','nfr',num_iint1*sizeofint) #endif do n=1,num_iint1f1 iint1(n) = iint1f1(n) nfr(n) = 1 enddo do n=1,num_iint1f2 iint1(n+num_iint1f1) = iint1f2(n) nfr(n+num_iint1f1) = 2 enddo ! ! Allocate stencil arrays with enough memory to hold both first and ! second layer fringe points ! ql1: quality using first-layer fringes ! ql2: quality using first- and second-layer fringes ! cd: cell difference ! d[xyz]t: stencil weights ! msh: donor mesh number ! ist: donor cell index ! ! if( num_iint1 .gt. 0 ) then allocate( ql1(num_iint1), stat=ierr) if( ierr /= 0 ) then write(lu_log,*) " ql1: allocation failed" call peg_error('bad allocate stat') endif #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','ql1',num_iint1*sizeofreal) #endif ql1 = -1.0 allocate( ql2(num_iint1), stat=ierr) if( ierr /= 0 ) then write(lu_log,*) " ql2: allocation failed" call peg_error('bad allocate stat') endif #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','ql2',num_iint1*sizeofreal) #endif ql2 = -1.0 allocate( cd(num_iint1), stat=ierr) if( ierr /= 0 ) then write(lu_log,*) " cd: allocation failed" call peg_error('bad allocate stat') endif #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','cd',num_iint1*sizeofreal) #endif cd = 1.0e30 allocate( dxt(num_iint1), stat=ierr) if( ierr /= 0 ) then write(lu_log,*) " dxt: allocation failed" call peg_error('bad allocate stat') endif #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','dxt',num_iint1*sizeofreal) #endif allocate( dyt(num_iint1), stat=ierr) if( ierr /= 0 ) then write(lu_log,*) " dyt: allocation failed" call peg_error('bad allocate stat') endif #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','dyt',num_iint1*sizeofreal) #endif allocate( dzt(num_iint1), stat=ierr) if( ierr /= 0 ) then write(lu_log,*) " dzt: allocation failed" call peg_error('bad allocate stat') endif #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','dzt',num_iint1*sizeofreal) #endif allocate( msh(num_iint1), stat=ierr) if( ierr /= 0 ) then write(lu_log,*) " msh: allocation failed" call peg_error('bad allocate stat') endif #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','msh',num_iint1*sizeofreal) #endif msh = 0 allocate( ist(num_iint1), stat=ierr) if( ierr /= 0 ) then write(lu_log,*) " ist: allocation failed" call peg_error('bad allocate stat') endif #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','ist',num_iint1*sizeofreal) #endif ! ! Loop through each possible donor mesh ! do m = 1,proc%num_secondary_1 m2 = mesh_no(proc%secondary_1(m)) ! ! Load stencils for this donor mesh ! call load_intrp( proc%primary_1, proc%secondary_1(m), & num_isten, isten ) ! call modify_interior( proc%secondary_1(m), & ! num_isten, isten ) if( num_isten .gt. 0 ) then allocate( iblank(mjmax(m2),mkmax(m2),mlmax(m2)), stat=ierr) if( ierr /= 0 ) then write(lu_log,*) " iblank: allocation failed" call peg_error('bad allocate stat') endif #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','iblank', & mjmax(m2)*mkmax(m2)*mlmax(m2)*sizeofint) #endif ! ! Load mesh coordinates for donor and recipient, hole points for ! donor, and lists of recipient points for donor and recipient meshes ! call load_x( proc%primary_1, x, y, z ) call load_x( proc%secondary_1(m), x2, y2, z2 ) call load_hole( proc%secondary_1(m), num_ihole, ihole ) ! call load_int1( proc%primary_1, num_iint1f1, iint1f1, & ! num_iint1f2, iint1f2 ) call load_int1( proc%secondary_1(m), num_iint1f1b, iint1f1b, & num_iint1f2b, iint1f2b ) ! ! Build an iblank array for this donor mesh ! 0 = hole point ! 1 = interior point ! -2 = first-layer fringe point ! -3 = second-layer fringe point ! iblank = 1 ! ! Composite hole ! do n = 1,num_ihole call itojkl( ihole(n), ji, ki, li, m2 ) iblank(ji,ki,li) = 0 enddo ! ! type='INTR' hole: treat the edges of the volume as single-fringe ! points (even though they are not explicitly added to the ! fringe-point list), and then blank out the interior ! do nr = 1,nregn if( rtype(nr) .eq. 'INTR' ) then if( ispartmr(nr) .eq. proc%secondary_1(m) ) then do nv = 1,nvol if( ispartr(nv) .eq. rname(nr) ) then jbeg = max(1, jrangev(1,nv)) jend = min(mjmax(m2), jrangev(2,nv)) kbeg = max(1, krangev(1,nv)) kend = min(mkmax(m2), krangev(2,nv)) lbeg = max(1, lrangev(1,nv)) lend = min(mlmax(m2), lrangev(2,nv)) do li=lbeg,lend do ki=kbeg,kend do ji=jbeg,jend iblank(ji,ki,li) = -2 enddo enddo enddo do li=lrangev(1,nv)+1,lrangev(2,nv)-1 do ki=krangev(1,nv)+1,krangev(2,nv)-1 do ji=jrangev(1,nv)+1,jrangev(2,nv)-1 iblank(ji,ki,li) = 0 enddo enddo enddo endif enddo endif endif enddo ! ! First layer fringes ! do n = 1,num_iint1f1b call itojkl( iint1f1b(n), ji, ki, li, m2 ) iblank(ji,ki,li) = -2 enddo ! ! Second layer fringes ! do n = 1,num_iint1f2b call itojkl( iint1f2b(n), ji, ki, li, m2 ) iblank(ji,ki,li) = -3 enddo ! ! If there is a projection for this mesh pair, load projection and ! shift the recipient mesh. ! if( project(m1,m2) ) then call load_project(proc%primary_1, proc%secondary_1(m), & num_project, i_project, & idir_project, x_project, & y_project, z_project ) if( num_project .gt. 0 ) then call shift_mesh(x,y,z,mjmax(m1),mkmax(m1),mlmax(m1), & m1,shift(m1,m2), & i_project, idir_project, & x_project,y_project,z_project,num_project) endif endif ! ! Compute the stencil weights, quality, and cell-difference ! parameter. Keep the stencil if it meets the quality criteria ! and if the cd is better than other donors ! call comp_level1( & m1, mjmax(m1), mkmax(m1), mlmax(m1), x, y, z, & m2, mjmax(m2), mkmax(m2), mlmax(m2), x2, y2, z2, & iblank, num_iint1, iint1, nfr, isten, & cd, ql1, ql2, dxt, dyt, dzt, msh, ist, & iblank_tmp ) ! ! un-shift mesh if projected. ! if( project(m1,m2) .and. num_project.gt.0 ) then call unshift_mesh(x,y,z,mjmax(m1),mkmax(m1),mlmax(m1)) endif deallocate(iblank, stat=ierr) if( ierr /= 0 ) call peg_error('bad deallocate stat for iblank') #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','iblank',-1) #endif endif enddo endif ! ! loop through each donor mesh and store the stencils ! do m = 1,proc%num_secondary_1 m2 = mesh_no(proc%secondary_1(m)) num_ilevel1 = 0 do n = 1,num_iint1 if( msh(n) .eq. m2 ) num_ilevel1 = num_ilevel1 + 1 enddo l1data%num_ilevel1 = num_ilevel1 call alloc_level1( l1data ) nl = 0 do n = 1,num_iint1 if( msh(n) .eq. m2 ) then nl = nl + 1 l1data%ibp(nl) = iint1(n) l1data%isp(nl) = ist(n) l1data%nfrl(nl) = nfr(n) l1data%dxint(nl) = dxt(n) l1data%dyint(nl) = dyt(n) l1data%dzint(nl) = dzt(n) l1data%qual1(nl) = ql1(n) l1data%qual2(nl) = ql2(n) l1data%cell_diff(nl) = cd(n) #ifdef DEBUG if ( n .eq. 123 ) then write(*,*) ' ' write(*,*) 'in mod_spec_level1: ' write(*,*) ' n, nl = ',n,nl write(*,*) ' primary mesh = ',proc%primary_1 write(*,*) ' secondary mesh = ',proc%secondary_1(m) write(*,*) ' intrp weights : ', & dxt(n),dyt(n),dzt(n) write(*,*) ' ql1,ql2,cd: ',ql1(n),ql2(n),cd(n) write(*,*) ' ist,msh: ',ist(n),msh(n) endif #endif endif enddo call store_level1( proc%primary_1, proc%secondary_1(m), l1data, 0 ) write(lu_log,1) num_ilevel1, proc%primary_1, proc%secondary_1(m) 1 format(' FOUND:',i10,' level1 interpolated points for:',/, & ' Recipient Mesh:',a40,/, & ' Donor Mesh:',a40,/) enddo ! ! find and store the orphan points ! num_iorphan = 0 do n = 1,num_iint1 if( msh(n) .eq. 0 ) num_iorphan = num_iorphan + 1 enddo orph_data%num_iorphan = num_iorphan call alloc_orphan( orph_data ) nl = 0 do n = 1,num_iint1 if( msh(n) .eq. 0 ) then nl = nl + 1 orph_data%iorphan(nl) = iint1(n) orph_data%nfrl(nl) = nfr(n) endif enddo call store_orphan( proc%primary_1, orph_data ) write(lu_log,2) num_iorphan,trim(proc%primary_1) 2 format(' FOUND:',i10,' ORPHAN points for:',/, & ' Mesh:',a40,/) ! ! Deallocate memory ! if( num_iint1 .gt. 0 ) then deallocate(ql1, stat=ierr) if( ierr /= 0 ) call peg_error('bad deallocate stat for ql1') #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','ql1',-1) #endif deallocate(ql2, stat=ierr) if( ierr /= 0 ) call peg_error('bad deallocate stat for ql2') #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','ql2',-1) #endif deallocate(cd, stat=ierr) if( ierr /= 0 ) call peg_error('bad deallocate stat for cd') #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','cd',-1) #endif deallocate(dxt, stat=ierr) if( ierr /= 0 ) call peg_error('bad deallocate stat for dxt') #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','dxt',-1) #endif deallocate(dyt, stat=ierr) if( ierr /= 0 ) call peg_error('bad deallocate stat for dyt') #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','dyt',-1) #endif deallocate(dzt, stat=ierr) if( ierr /= 0 ) call peg_error('bad deallocate stat for dzt') #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','dzt',-1) #endif deallocate(msh, stat=ierr) if( ierr /= 0 ) call peg_error('bad deallocate stat for msh') #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','msh',-1) #endif deallocate(ist, stat=ierr) if( ierr /= 0 ) call peg_error('bad deallocate stat for ist') #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','ist',-1) #endif endif deallocate(iblank_tmp, iint1, nfr, stat=ierr) if( ierr /= 0 ) call peg_error('bad deallocate stat for iblank_tmp') #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','iblank_tmp',-1) call mem_track('execute_spec_level1','iint1',-1) call mem_track('execute_spec_level1','nfr',-1) #endif deallocate( l1data%ibp, l1data%isp, l1data%nfrl, & l1data%dxint, l1data%dyint, l1data%dzint, & l1data%qual1, l1data%qual2, l1data%cell_diff, stat=ierr) if( ierr /= 0 ) call peg_error('bad deallocate stat for l1data') #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','ibp',-1) call mem_track('execute_spec_level1','isp',-1) call mem_track('execute_spec_level1','nfrl',-1) call mem_track('execute_spec_level1','dxint',-1) call mem_track('execute_spec_level1','dyint',-1) call mem_track('execute_spec_level1','dzint',-1) call mem_track('execute_spec_level1','qual1',-1) call mem_track('execute_spec_level1','qual2',-1) call mem_track('execute_spec_level1','cell_diff',-1) #endif deallocate(orph_data%iorphan, orph_data%nfrl, stat=ierr) if( ierr /= 0 ) call peg_error('bad deallocate stat for orph_data') #ifdef MEMORY_TRACKING call mem_track('execute_spec_level1','iorphan',-1) call mem_track('execute_spec_level1','nfrl',-1) #endif return end subroutine execute_spec_level1 end module mod_spec_level1