[cig-commits] r20561 - in seismo/2D/SPECFEM2D/trunk/src: meshfem2D specfem2D

xie.zhinan at geodynamics.org xie.zhinan at geodynamics.org
Wed Aug 8 19:32:21 PDT 2012


Author: xie.zhinan
Date: 2012-08-08 19:32:20 -0700 (Wed, 08 Aug 2012)
New Revision: 20561

Modified:
   seismo/2D/SPECFEM2D/trunk/src/meshfem2D/read_parameter_file.F90
   seismo/2D/SPECFEM2D/trunk/src/meshfem2D/save_databases.f90
   seismo/2D/SPECFEM2D/trunk/src/specfem2D/pml_init.F90
   seismo/2D/SPECFEM2D/trunk/src/specfem2D/read_databases.f90
   seismo/2D/SPECFEM2D/trunk/src/specfem2D/specfem2D.F90
Log:
add first edition of supporting simulation with CPML when read_external_mesh is true


Modified: seismo/2D/SPECFEM2D/trunk/src/meshfem2D/read_parameter_file.F90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/src/meshfem2D/read_parameter_file.F90	2012-08-09 01:22:39 UTC (rev 20560)
+++ seismo/2D/SPECFEM2D/trunk/src/meshfem2D/read_parameter_file.F90	2012-08-09 02:32:20 UTC (rev 20561)
@@ -58,7 +58,8 @@
   logical :: SAVE_FORWARD,read_external_mesh
 
   character(len=256) :: mesh_file, nodes_coords_file, materials_file, &
-                        free_surface_file, absorbing_surface_file
+                        free_surface_file, absorbing_surface_file,&
+                        CPML_element_file
   character(len=256)  :: tangential_detection_curve_file
 
   ! variables used for partitioning
@@ -499,6 +500,9 @@
   call read_value_string_p(absorbing_surface_file, 'mesher.absorbing_surface_file')
   if(err_occurred() /= 0) stop 'error reading parameter 56 in Par_file'
 
+  call read_value_string_p(CPML_element_file, 'mesher.CPML_element_file')
+  if(err_occurred() /= 0) stop 'error reading parameter 56 in Par_file'
+
   call read_value_string_p(tangential_detection_curve_file, 'mesher.tangential_detection_curve_file')
   if(err_occurred() /= 0) stop 'error reading parameter 57 in Par_file'
 

Modified: seismo/2D/SPECFEM2D/trunk/src/meshfem2D/save_databases.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/src/meshfem2D/save_databases.f90	2012-08-09 01:22:39 UTC (rev 20560)
+++ seismo/2D/SPECFEM2D/trunk/src/meshfem2D/save_databases.f90	2012-08-09 02:32:20 UTC (rev 20561)
@@ -118,6 +118,14 @@
     write(15,*) 'PML_BOUNDARY_CONDITIONS'
     write(15,*) PML_BOUNDARY_CONDITIONS
 
+    write(15,*) 'read_external_mesh'
+    write(15,*) read_external_mesh
+
+    if(read_external_mesh) then
+    write(15,*) 'CPML_element_file'
+    write(15,"(a256)") trim(CPML_element_file)
+    endif
+
     write(15,*) 'NELEM_PML_THICKNESS'
     write(15,*) NELEM_PML_THICKNESS
 

Modified: seismo/2D/SPECFEM2D/trunk/src/specfem2D/pml_init.F90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/src/specfem2D/pml_init.F90	2012-08-09 01:22:39 UTC (rev 20560)
+++ seismo/2D/SPECFEM2D/trunk/src/specfem2D/pml_init.F90	2012-08-09 02:32:20 UTC (rev 20561)
@@ -43,7 +43,8 @@
 
   subroutine pml_init(nspec,nglob,anyabs,ibool,nelemabs,codeabs,numabs,&
                     nspec_PML,is_PML,which_PML_elem,spec_to_PML, &
-                    icorner_iglob,NELEM_PML_THICKNESS)
+                    icorner_iglob,NELEM_PML_THICKNESS,&
+                  read_external_mesh,CPML_element_file)
 
 
   implicit none
@@ -67,8 +68,16 @@
   integer, dimension(nspec) :: spec_to_PML
   logical, dimension(nspec) :: is_PML
 
+!! DK DK for CPML_element_file
+  logical :: read_external_mesh          
+  character(len=256)  :: CPML_element_file 
+  integer :: ier,ispec_CPML
+  integer :: ispec_temp
+  integer, dimension(:), allocatable :: region_CPML
+
   !!!detection of PML elements
 
+  if(.not. read_external_mesh)then
   nspec_PML = 0
 
      !ibound is the side we are looking (bottom, right, top or left)
@@ -162,6 +171,85 @@
 
      write(IOUT,*) "number of PML spectral elements :", nspec_PML
 
+     endif
+
+!!!!!!when read_external_mesh read the element in PML
+  if(read_external_mesh)then
+  is_PML(:) = .false.
+  which_PML_elem(:,:) = .false.
+
+  open(unit=9999, file=trim(CPML_element_file), form='formatted' , status='old', action='read',iostat=ier)
+  if( ier /= 0 ) then
+    print*,'error opening file: ',trim(CPML_element_file)
+    stop 'error read external CPML element file'
+  endif  
+
+  read(9999,*)nspec_PML 
+  allocate(region_CPML(nspec_PML))
+
+  do ispec_CPML=1,nspec_PML 
+     read(9999,*)ispec_temp,region_CPML(ispec_CPML)
+     is_PML(ispec_temp)=.true.
+     spec_to_PML(ispec_temp)=ispec_CPML
+  enddo
+
+  do ispec=1,nspec
+     if(is_PML(ispec))then
+     ispec_CPML=spec_to_PML(ispec)
+       if(region_CPML(ispec_CPML)==1)then
+         which_PML_elem(ILEFT,ispec)   = .true.
+         which_PML_elem(IRIGHT,ispec)  = .false.
+         which_PML_elem(ITOP,ispec)    = .false.
+         which_PML_elem(IBOTTOM,ispec) = .false.
+       elseif(region_CPML(ispec_CPML)==2)then
+         which_PML_elem(ILEFT,ispec)   = .false.
+         which_PML_elem(IRIGHT,ispec)  = .true.
+         which_PML_elem(ITOP,ispec)    = .false.
+         which_PML_elem(IBOTTOM,ispec) = .false.
+       elseif(region_CPML(ispec_CPML)==4)then
+         which_PML_elem(ILEFT,ispec)   = .false.
+         which_PML_elem(IRIGHT,ispec)  = .false.
+         which_PML_elem(ITOP,ispec)    = .true.
+         which_PML_elem(IBOTTOM,ispec) = .false.
+       elseif(region_CPML(ispec_CPML)==5)then
+         which_PML_elem(ILEFT,ispec)   = .true.
+         which_PML_elem(IRIGHT,ispec)  = .false.
+         which_PML_elem(ITOP,ispec)    = .true.
+         which_PML_elem(IBOTTOM,ispec) = .false.
+       elseif(region_CPML(ispec_CPML)==6)then
+         which_PML_elem(ILEFT,ispec)   = .false.
+         which_PML_elem(IRIGHT,ispec)  = .true.
+         which_PML_elem(ITOP,ispec)    = .true.
+         which_PML_elem(IBOTTOM,ispec) = .false.
+       elseif(region_CPML(ispec_CPML)==8)then
+         which_PML_elem(ILEFT,ispec)   = .false.
+         which_PML_elem(IRIGHT,ispec)  = .false.
+         which_PML_elem(ITOP,ispec)    = .false.
+         which_PML_elem(IBOTTOM,ispec) = .true.
+       elseif(region_CPML(ispec_CPML)==9)then
+         which_PML_elem(ILEFT,ispec)   = .true.
+         which_PML_elem(IRIGHT,ispec)  = .false.
+         which_PML_elem(ITOP,ispec)    = .false.
+         which_PML_elem(IBOTTOM,ispec) = .true.
+       elseif(region_CPML(ispec_CPML)==10)then
+         which_PML_elem(ILEFT,ispec)   = .false.
+         which_PML_elem(IRIGHT,ispec)  = .true.
+         which_PML_elem(ITOP,ispec)    = .false.
+         which_PML_elem(IBOTTOM,ispec) = .true.
+       else
+         which_PML_elem(ILEFT,ispec)   = .false.
+         which_PML_elem(IRIGHT,ispec)  = .false.
+         which_PML_elem(ITOP,ispec)    = .false.
+         which_PML_elem(IBOTTOM,ispec) = .false.
+       endif
+     endif
+  enddo
+
+     write(IOUT,*) "number of PML spectral elements :", nspec_PML
+
+  endif
+!!!!!!when read_external_mesh read the element in PML
+
   end subroutine pml_init
 
 !

Modified: seismo/2D/SPECFEM2D/trunk/src/specfem2D/read_databases.f90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/src/specfem2D/read_databases.f90	2012-08-09 01:22:39 UTC (rev 20560)
+++ seismo/2D/SPECFEM2D/trunk/src/specfem2D/read_databases.f90	2012-08-09 02:32:20 UTC (rev 20561)
@@ -58,7 +58,8 @@
                   DRAW_SOURCES_AND_RECEIVERS,Q0,freq0,p_sv,NSTEP,deltat,NSOURCES, &
                   factor_subsample_image,USE_SNAPSHOT_NUMBER_IN_FILENAME,DRAW_WATER_IN_BLUE,US_LETTER, &
                   POWER_DISPLAY_COLOR,PERFORM_CUTHILL_MCKEE,SU_FORMAT,USER_T0,time_stepping_scheme,&
-                  ADD_SPRING_TO_STACEY,ADD_PERIODIC_CONDITIONS,PERIODIC_horiz_dist,PERIODIC_DETECT_TOL)
+                  ADD_SPRING_TO_STACEY,ADD_PERIODIC_CONDITIONS,PERIODIC_horiz_dist,PERIODIC_DETECT_TOL,&
+                  read_external_mesh,CPML_element_file)
 
 ! starts reading in parameters from input Database file
 
@@ -130,6 +131,10 @@
 !! DK DK grid point detection tolerance for periodic conditions
   double precision :: PERIODIC_DETECT_TOL
 
+!! DK DK for CPML_element_file
+  logical :: read_external_mesh          
+  character(len=256)  :: CPML_element_file  
+
   ! local parameters
   integer :: ier
   character(len=80) :: datlin
@@ -183,6 +188,14 @@
   read(IIN,"(a80)") datlin
   read(IIN,*) PML_BOUNDARY_CONDITIONS
 
+  read(IIN,"(a80)") datlin    
+  read(IIN,*) read_external_mesh  
+
+  if(read_external_mesh)then     
+  read(IIN,"(a80)") datlin      
+  read(IIN,"(a256)") CPML_element_file  
+  endif                          
+
   read(IIN,"(a80)") datlin
   read(IIN,*) NELEM_PML_THICKNESS
 

Modified: seismo/2D/SPECFEM2D/trunk/src/specfem2D/specfem2D.F90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/src/specfem2D/specfem2D.F90	2012-08-09 01:22:39 UTC (rev 20560)
+++ seismo/2D/SPECFEM2D/trunk/src/specfem2D/specfem2D.F90	2012-08-09 02:32:20 UTC (rev 20561)
@@ -561,6 +561,10 @@
     plot_lowerleft_corner_only,add_Bielak_conditions,output_energy,READ_EXTERNAL_SEP_FILE, &
     output_wavefield_dumps,use_binary_for_wavefield_dumps,PML_BOUNDARY_CONDITIONS
 
+!! DK DK for CPML_element_file
+  logical :: read_external_mesh           
+  character(len=256)  :: CPML_element_file  
+
   double precision :: cutsnaps,sizemax_arrows,anglerec,xirec,gammarec
 
 ! for absorbing and acoustic free surface conditions
@@ -1041,7 +1045,8 @@
                   Q0,freq0,p_sv,NSTEP,deltat,NSOURCES, &
                   factor_subsample_image,USE_SNAPSHOT_NUMBER_IN_FILENAME,DRAW_WATER_IN_BLUE,US_LETTER, &
                   POWER_DISPLAY_COLOR,PERFORM_CUTHILL_MCKEE,SU_FORMAT,USER_T0, time_stepping_scheme, &
-                  ADD_SPRING_TO_STACEY,ADD_PERIODIC_CONDITIONS,PERIODIC_horiz_dist,PERIODIC_DETECT_TOL)
+                  ADD_SPRING_TO_STACEY,ADD_PERIODIC_CONDITIONS,PERIODIC_horiz_dist,PERIODIC_DETECT_TOL,&
+                  read_external_mesh,CPML_element_file) 
   if(nproc_read_from_database < 1) stop 'should have nproc_read_from_database >= 1'
   if(SIMULATION_TYPE == 2 .and.(time_stepping_scheme == 2 .or. time_stepping_scheme == 3)) &
                                   stop 'RK and LDDRK time scheme not supported for adjoint inversion'
@@ -1081,7 +1086,8 @@
                       Q0,freq0,p_sv,NSTEP,deltat,NSOURCES, &
                       factor_subsample_image,USE_SNAPSHOT_NUMBER_IN_FILENAME,DRAW_WATER_IN_BLUE,US_LETTER, &
                       POWER_DISPLAY_COLOR,PERFORM_CUTHILL_MCKEE,SU_FORMAT,USER_T0, time_stepping_scheme, &
-                      ADD_SPRING_TO_STACEY,ADD_PERIODIC_CONDITIONS,PERIODIC_horiz_dist,PERIODIC_DETECT_TOL)
+                      ADD_SPRING_TO_STACEY,ADD_PERIODIC_CONDITIONS,PERIODIC_horiz_dist,PERIODIC_DETECT_TOL,&
+                      read_external_mesh,CPML_element_file)
 
   !
   !--- source information
@@ -2820,9 +2826,14 @@
       is_PML(:) = .false.
       which_PML_elem(:,:) = .false.
 
+
+      write(*,*)read_external_mesh,'read_external_mesh'
+      write(*,*)CPML_element_file,'CPML_element_file'
+
       call pml_init(nspec,nglob,anyabs,ibool,nelemabs,codeabs,numabs,&
                   nspec_PML,is_PML,which_PML_elem,spec_to_PML, &
-                  icorner_iglob,NELEM_PML_THICKNESS)
+                  icorner_iglob,NELEM_PML_THICKNESS,&
+                  read_external_mesh,CPML_element_file)
 
       deallocate(icorner_iglob)
 



More information about the CIG-COMMITS mailing list