diff --git a/src/specfem2D/define_external_model_from_tomo_file.f90 b/src/specfem2D/define_external_model_from_tomo_file.f90 index 6a8e91d3d..19d3aaff9 100644 --- a/src/specfem2D/define_external_model_from_tomo_file.f90 +++ b/src/specfem2D/define_external_model_from_tomo_file.f90 @@ -367,7 +367,7 @@ subroutine define_external_model_from_tomo_file(rhoext,vpext,vsext, & c33ext,c35ext,c55ext ! local parameters - integer :: i,j,ispec,iglob + integer :: i,j,ispec,iglob,ival double precision :: xmesh,zmesh double precision :: rho_final double precision :: vp_final,vs_final @@ -407,9 +407,10 @@ subroutine define_external_model_from_tomo_file(rhoext,vpext,vsext, & do ispec = 1,nspec do j = 1,NGLLZ do i = 1,NGLLX - iglob = ibool(i,j,ispec) + ! determine material properties if (kmato(ispec) == tomo_material) then ! If the material has been set to < 0 on the Par_file + iglob = ibool(i,j,ispec) xmesh = coord(1,iglob) zmesh = coord(2,iglob) @@ -512,10 +513,10 @@ subroutine define_external_model_from_tomo_file(rhoext,vpext,vsext, & call synchronize_all() ! collect totals - iglob = npoint_tomo - call max_all_i(iglob,npoint_tomo) - iglob = npoint_internal - call max_all_i(iglob,npoint_internal) + ival = npoint_tomo + call max_all_i(ival,npoint_tomo) + ival = npoint_internal + call max_all_i(ival,npoint_internal) ! user output if (myrank == 0) then diff --git a/src/specfem2D/setup_mesh.F90 b/src/specfem2D/setup_mesh.F90 index 2cfd3a505..25f47c238 100644 --- a/src/specfem2D/setup_mesh.F90 +++ b/src/specfem2D/setup_mesh.F90 @@ -561,6 +561,8 @@ subroutine setup_mesh_material_properties() double precision :: c33min_glob,c33max_glob,c35min_glob,c35max_glob double precision :: c55min_glob,c55max_glob + logical :: has_poroelasticity, has_anisotropy + ! vtk output character(len=MAX_STRING_LEN) :: filename,prname double precision,dimension(:,:,:),allocatable :: tmp_store @@ -945,16 +947,25 @@ subroutine setup_mesh_material_properties() tmp_val = maxval(qmu_attenuation_store) call max_all_dp(tmp_val, qmumax_glob) - if (any_poroelastic) then + ! poroelasticity + ! check if any poroelastic elements in domains + call any_all_l(any_poroelastic,has_poroelasticity) + if (has_poroelasticity) then ! vpII - call min_all_dp(vpIImin_glob, vpIImin_glob) - call max_all_dp(vpIImax_glob, vpIImax_glob) + tmp_val = vpIImin_glob + call min_all_dp(tmp_val, vpIImin_glob) + tmp_val = vpIImax_glob + call max_all_dp(tmp_val, vpIImax_glob) ! phi - call min_all_dp(phimin_glob, phimin_glob) - call max_all_dp(phimax_glob, phimax_glob) + tmp_val = phimin_glob + call min_all_dp(tmp_val, phimin_glob) + tmp_val = phimax_glob + call max_all_dp(tmp_val, phimax_glob) endif - if (any_anisotropy) then + ! anisotropy + call any_all_l(any_anisotropy,has_anisotropy) + if (has_anisotropy) then ! c11 tmp_val = minval(c11store) call min_all_dp(tmp_val, c11min_glob) @@ -1017,13 +1028,13 @@ subroutine setup_mesh_material_properties() write(IMAIN,*) ' vs : min/max = ',sngl(vsmin_glob),'/',sngl(vsmax_glob) write(IMAIN,*) ' vp : min/max = ',sngl(vpmin_glob),'/',sngl(vpmax_glob) write(IMAIN,*) - if (any_poroelastic) then + if (has_poroelasticity) then write(IMAIN,*) ' poroelasticity:' write(IMAIN,*) ' vpII: min/max = ',sngl(vpIImin_glob),'/',sngl(vpIImax_glob) write(IMAIN,*) ' phi : min/max = ',sngl(phimin_glob),'/',sngl(phimax_glob) write(IMAIN,*) endif - if (any_anisotropy) then + if (has_anisotropy) then write(IMAIN,*) ' anisotropy:' write(IMAIN,*) ' c11: min/max = ',sngl(c11min_glob),'/',sngl(c11max_glob) write(IMAIN,*) ' c12: min/max = ',sngl(c12min_glob),'/',sngl(c12max_glob) diff --git a/src/specfem2D/specfem2D_par.f90 b/src/specfem2D/specfem2D_par.f90 index 491f8797a..00cab2d0a 100644 --- a/src/specfem2D/specfem2D_par.f90 +++ b/src/specfem2D/specfem2D_par.f90 @@ -423,8 +423,11 @@ module specfem_par ! anisotropy integer :: nspec_aniso + ! local flag if any anisotropic element is in this slice logical :: any_anisotropy + ! local flag if all elements are anisotropic in this slice logical :: all_anisotropic + logical, dimension(:), allocatable :: ispec_is_anisotropic ! inverse mass matrices