[cig-commits] r19799 - seismo/3D/CPML/trunk
dkomati1 at geodynamics.org
dkomati1 at geodynamics.org
Sun Mar 18 11:22:12 PDT 2012
Author: dkomati1
Date: 2012-03-18 11:22:12 -0700 (Sun, 18 Mar 2012)
New Revision: 19799
Modified:
seismo/3D/CPML/trunk/seismic_ADEPML_2D_RK4_eighth_order.f90
seismo/3D/CPML/trunk/seismic_CPML_2D_anisotropic.f90
seismo/3D/CPML/trunk/seismic_CPML_2D_isotropic_fourth_order.f90
seismo/3D/CPML/trunk/seismic_CPML_2D_isotropic_second_order.f90
seismo/3D/CPML/trunk/seismic_CPML_2D_poroelastic_fourth_order.f90
seismo/3D/CPML/trunk/seismic_CPML_3D_isotropic_MPI_OpenMP.f90
seismo/3D/CPML/trunk/seismic_CPML_3D_viscoelastic_MPI.f90
Log:
made sure alpha() is never zero to avoid having an infinite static gain in the Butterworth filter
Modified: seismo/3D/CPML/trunk/seismic_ADEPML_2D_RK4_eighth_order.f90
===================================================================
--- seismo/3D/CPML/trunk/seismic_ADEPML_2D_RK4_eighth_order.f90 2012-03-18 16:41:52 UTC (rev 19798)
+++ seismo/3D/CPML/trunk/seismic_ADEPML_2D_RK4_eighth_order.f90 2012-03-18 18:22:12 UTC (rev 19799)
@@ -286,16 +286,11 @@
! Kappa must be strong enough to absorb energy and low enough to avoid
! reflections.
-! Alpha1 must be low to absron energy and high enough to have efficiency on
+! Alpha1 must be low to absorb energy and high enough to have efficiency on
! grazing incident waves.
-! double precision, parameter :: K_MAX_PML_1 = 15.d0
double precision, parameter :: K_MAX_PML_1 = 7.d0
double precision, parameter :: ALPHA_MAX_PML_1 = 2.d0*PI*(f0/2.d0)
-! double precision, parameter :: K_MAX_PML_2 = K_MAX_PML_1 / 15.d0
- double precision, parameter :: K_MAX_PML_2 = K_MAX_PML_1
- double precision, parameter :: ALPHA_MAX_PML_2 = ALPHA_MAX_PML_1 / 5.d0
-
! arrays for the memory variables
! could declare these arrays in PML only to save a lot of memory, but proof of concept only here
! We have as many memory variables as the number of frequency shift poles in the CPML
@@ -507,7 +502,7 @@
d_x_1(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x_1(i) = 1.d0 + (K_MAX_PML_1 - 1.d0) * abscissa_normalized**NPOWER2
- alpha_x_1(i) = ALPHA_MAX_PML_1 * (1.d0 - abscissa_normalized)
+ alpha_x_1(i) = ALPHA_MAX_PML_1 * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML_1
endif
! define damping profile at half the grid points
@@ -517,7 +512,7 @@
d_x_half_1(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x_half_1(i) = 1.d0 + (K_MAX_PML_1 - 1.d0) * abscissa_normalized**NPOWER2
- alpha_x_half_1(i) = ALPHA_MAX_PML_1 * (1.d0 - abscissa_normalized)
+ alpha_x_half_1(i) = ALPHA_MAX_PML_1 * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML_1
endif
endif
@@ -532,7 +527,7 @@
d_x_1(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x_1(i) = 1.d0 + (K_MAX_PML_1 - 1.d0) * abscissa_normalized**NPOWER2
- alpha_x_1(i) = ALPHA_MAX_PML_1 * (1.d0 - abscissa_normalized)
+ alpha_x_1(i) = ALPHA_MAX_PML_1 * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML_1
endif
! define damping profile at half the grid points
@@ -542,7 +537,7 @@
d_x_half_1(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x_half_1(i) = 1.d0 + (K_MAX_PML_1 - 1.d0) * abscissa_normalized**NPOWER2
- alpha_x_half_1(i) = ALPHA_MAX_PML_1 * (1.d0 - abscissa_normalized)
+ alpha_x_half_1(i) = ALPHA_MAX_PML_1 * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML_1
endif
endif
@@ -601,7 +596,7 @@
d_y_1(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y_1(j) = 1.d0 + (K_MAX_PML_1 - 1.d0) * abscissa_normalized**NPOWER2
- alpha_y_1(j) = ALPHA_MAX_PML_1 * (1.d0 - abscissa_normalized)
+ alpha_y_1(j) = ALPHA_MAX_PML_1 * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML_1
endif
! define damping profile at half the grid points
@@ -611,7 +606,7 @@
d_y_half_1(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y_half_1(j) = 1.d0 + (K_MAX_PML_1 - 1.d0) * abscissa_normalized**NPOWER2
- alpha_y_half_1(j) = ALPHA_MAX_PML_1 * (1.d0 - abscissa_normalized)
+ alpha_y_half_1(j) = ALPHA_MAX_PML_1 * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML_1
endif
endif
@@ -626,7 +621,7 @@
d_y_1(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y_1(j) = 1.d0 + (K_MAX_PML_1 - 1.d0) * abscissa_normalized**NPOWER2
- alpha_y_1(j) = ALPHA_MAX_PML_1 * (1.d0 - abscissa_normalized)
+ alpha_y_1(j) = ALPHA_MAX_PML_1 * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML_1
endif
! define damping profile at half the grid points
@@ -636,7 +631,7 @@
d_y_half_1(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y_half_1(j) = 1.d0 + (K_MAX_PML_1 - 1.d0) * abscissa_normalized**NPOWER2
- alpha_y_half_1(j) = ALPHA_MAX_PML_1 * (1.d0 - abscissa_normalized)
+ alpha_y_half_1(j) = ALPHA_MAX_PML_1 * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML_1
endif
endif
Modified: seismo/3D/CPML/trunk/seismic_CPML_2D_anisotropic.f90
===================================================================
--- seismo/3D/CPML/trunk/seismic_CPML_2D_anisotropic.f90 2012-03-18 16:41:52 UTC (rev 19798)
+++ seismo/3D/CPML/trunk/seismic_CPML_2D_anisotropic.f90 2012-03-18 18:22:12 UTC (rev 19799)
@@ -436,7 +436,7 @@
d_x(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -446,7 +446,7 @@
d_x_half(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x_half(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x_half(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x_half(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
@@ -461,7 +461,7 @@
d_x(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -471,7 +471,7 @@
d_x_half(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x_half(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x_half(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x_half(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
@@ -511,7 +511,7 @@
d_y(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -521,7 +521,7 @@
d_y_half(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y_half(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y_half(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y_half(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
@@ -536,7 +536,7 @@
d_y(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -546,7 +546,7 @@
d_y_half(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y_half(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y_half(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y_half(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
Modified: seismo/3D/CPML/trunk/seismic_CPML_2D_isotropic_fourth_order.f90
===================================================================
--- seismo/3D/CPML/trunk/seismic_CPML_2D_isotropic_fourth_order.f90 2012-03-18 16:41:52 UTC (rev 19798)
+++ seismo/3D/CPML/trunk/seismic_CPML_2D_isotropic_fourth_order.f90 2012-03-18 18:22:12 UTC (rev 19799)
@@ -349,7 +349,7 @@
d_x(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -359,7 +359,7 @@
d_x_half(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x_half(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x_half(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x_half(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
@@ -374,7 +374,7 @@
d_x(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -384,7 +384,7 @@
d_x_half(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x_half(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x_half(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x_half(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
@@ -424,7 +424,7 @@
d_y(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -434,7 +434,7 @@
d_y_half(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y_half(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y_half(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y_half(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
@@ -449,7 +449,7 @@
d_y(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -459,7 +459,7 @@
d_y_half(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y_half(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y_half(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y_half(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
Modified: seismo/3D/CPML/trunk/seismic_CPML_2D_isotropic_second_order.f90
===================================================================
--- seismo/3D/CPML/trunk/seismic_CPML_2D_isotropic_second_order.f90 2012-03-18 16:41:52 UTC (rev 19798)
+++ seismo/3D/CPML/trunk/seismic_CPML_2D_isotropic_second_order.f90 2012-03-18 18:22:12 UTC (rev 19799)
@@ -348,7 +348,7 @@
d_x(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -358,7 +358,7 @@
d_x_half(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x_half(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x_half(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x_half(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
@@ -373,7 +373,7 @@
d_x(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -383,7 +383,7 @@
d_x_half(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x_half(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x_half(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x_half(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
@@ -423,7 +423,7 @@
d_y(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -433,7 +433,7 @@
d_y_half(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y_half(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y_half(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y_half(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
@@ -448,7 +448,7 @@
d_y(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -458,7 +458,7 @@
d_y_half(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y_half(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y_half(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y_half(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
Modified: seismo/3D/CPML/trunk/seismic_CPML_2D_poroelastic_fourth_order.f90
===================================================================
--- seismo/3D/CPML/trunk/seismic_CPML_2D_poroelastic_fourth_order.f90 2012-03-18 16:41:52 UTC (rev 19798)
+++ seismo/3D/CPML/trunk/seismic_CPML_2D_poroelastic_fourth_order.f90 2012-03-18 18:22:12 UTC (rev 19799)
@@ -293,9 +293,7 @@
! power to compute d0 profile
double precision, parameter :: NPOWER = 2.d0
-! double precision, parameter :: K_MAX_PML = 7.d0 ! from Gedney page 8.11
double precision, parameter :: K_MAX_PML = 1.d0 ! from Gedney page 8.11
-! double precision, parameter :: ALPHA_MAX_PML = 0.05d0 ! from Gedney page 8.22
double precision, parameter :: ALPHA_MAX_PML = 2.d0*PI*(f0/2.d0) ! from Festa and Vilotte
! 2D arrays for the memory variables
@@ -456,7 +454,7 @@
d_x(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -466,7 +464,7 @@
d_x_half_x(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x_half_x(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x_half_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x_half_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
@@ -481,7 +479,7 @@
d_x(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -491,7 +489,7 @@
d_x_half_x(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x_half_x(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x_half_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x_half_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
@@ -532,7 +530,7 @@
d_y(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -542,7 +540,7 @@
d_y_half_y(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y_half_y(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y_half_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y_half_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
@@ -557,7 +555,7 @@
d_y(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -567,7 +565,7 @@
d_y_half_y(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y_half_y(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y_half_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y_half_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
Modified: seismo/3D/CPML/trunk/seismic_CPML_3D_isotropic_MPI_OpenMP.f90
===================================================================
--- seismo/3D/CPML/trunk/seismic_CPML_3D_isotropic_MPI_OpenMP.f90 2012-03-18 16:41:52 UTC (rev 19798)
+++ seismo/3D/CPML/trunk/seismic_CPML_3D_isotropic_MPI_OpenMP.f90 2012-03-18 18:22:12 UTC (rev 19799)
@@ -467,7 +467,7 @@
d_x(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -477,7 +477,7 @@
d_x_half(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x_half(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x_half(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x_half(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
@@ -492,7 +492,7 @@
d_x(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -502,7 +502,7 @@
d_x_half(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x_half(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x_half(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x_half(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
@@ -542,7 +542,7 @@
d_y(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -552,7 +552,7 @@
d_y_half(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y_half(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y_half(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y_half(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
@@ -567,7 +567,7 @@
d_y(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -577,7 +577,7 @@
d_y_half(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y_half(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y_half(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y_half(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
@@ -613,7 +613,7 @@
d_z(k) = d0_z * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_z(k) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_z(k) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_z(k) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -623,7 +623,7 @@
d_z_half(k) = d0_z * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_z_half(k) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_z_half(k) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_z_half(k) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
@@ -638,7 +638,7 @@
d_z(k) = d0_z * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_z(k) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_z(k) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_z(k) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -648,7 +648,7 @@
d_z_half(k) = d0_z * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_z_half(k) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_z_half(k) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_z_half(k) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
Modified: seismo/3D/CPML/trunk/seismic_CPML_3D_viscoelastic_MPI.f90
===================================================================
--- seismo/3D/CPML/trunk/seismic_CPML_3D_viscoelastic_MPI.f90 2012-03-18 16:41:52 UTC (rev 19798)
+++ seismo/3D/CPML/trunk/seismic_CPML_3D_viscoelastic_MPI.f90 2012-03-18 18:22:12 UTC (rev 19799)
@@ -216,7 +216,6 @@
double precision, parameter :: NPOWER = 2.d0
double precision, parameter :: K_MAX_PML = 7.d0 ! from Gedney page 8.11
-! double precision, parameter :: ALPHA_MAX_PML = 0.d0 ! from Festa and Vilotte
double precision, parameter :: ALPHA_MAX_PML = 2.d0*PI*(f0/2.d0) ! from Festa and Vilotte
! arrays for the memory variables
@@ -547,7 +546,7 @@
d_x(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -557,7 +556,7 @@
d_x_half(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x_half(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x_half(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x_half(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
@@ -572,7 +571,7 @@
d_x(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -582,7 +581,7 @@
d_x_half(i) = d0_x * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_x_half(i) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_x_half(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_x_half(i) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
@@ -622,7 +621,7 @@
d_y(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -632,7 +631,7 @@
d_y_half(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y_half(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y_half(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y_half(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
@@ -647,7 +646,7 @@
d_y(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -657,7 +656,7 @@
d_y_half(j) = d0_y * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_y_half(j) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_y_half(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_y_half(j) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
@@ -693,7 +692,7 @@
d_z(k) = d0_z * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_z(k) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_z(k) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_z(k) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -703,7 +702,7 @@
d_z_half(k) = d0_z * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_z_half(k) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_z_half(k) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_z_half(k) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
@@ -718,7 +717,7 @@
d_z(k) = d0_z * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_z(k) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_z(k) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_z(k) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
! define damping profile at half the grid points
@@ -728,7 +727,7 @@
d_z_half(k) = d0_z * abscissa_normalized**NPOWER
! this taken from Gedney page 8.2
K_z_half(k) = 1.d0 + (K_MAX_PML - 1.d0) * abscissa_normalized**NPOWER
- alpha_z_half(k) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized)
+ alpha_z_half(k) = ALPHA_MAX_PML * (1.d0 - abscissa_normalized) + 0.1d0 * ALPHA_MAX_PML
endif
endif
More information about the CIG-COMMITS
mailing list