[cig-commits] r18323 - in seismo/3D/FAULT_SOURCE/branches/new_fault_db: decompose_mesh_SCOTCH src
ampuero at geodynamics.org
ampuero at geodynamics.org
Thu May 5 17:24:52 PDT 2011
Author: ampuero
Date: 2011-05-05 17:24:51 -0700 (Thu, 05 May 2011)
New Revision: 18323
Modified:
seismo/3D/FAULT_SOURCE/branches/new_fault_db/decompose_mesh_SCOTCH/decompose_mesh_SCOTCH.f90
seismo/3D/FAULT_SOURCE/branches/new_fault_db/decompose_mesh_SCOTCH/fault_scotch.f90
seismo/3D/FAULT_SOURCE/branches/new_fault_db/src/create_regions_mesh.f90
seismo/3D/FAULT_SOURCE/branches/new_fault_db/src/fault_object.f90
Log:
fixed and cleaned up scotch and src
Modified: seismo/3D/FAULT_SOURCE/branches/new_fault_db/decompose_mesh_SCOTCH/decompose_mesh_SCOTCH.f90
===================================================================
--- seismo/3D/FAULT_SOURCE/branches/new_fault_db/decompose_mesh_SCOTCH/decompose_mesh_SCOTCH.f90 2011-05-05 19:52:23 UTC (rev 18322)
+++ seismo/3D/FAULT_SOURCE/branches/new_fault_db/decompose_mesh_SCOTCH/decompose_mesh_SCOTCH.f90 2011-05-06 00:24:51 UTC (rev 18323)
@@ -90,7 +90,7 @@
!----------------------------------------------------------------------------------------------
subroutine read_mesh_files
- ! sets number of nodes per element (Percy : esize = 8 defined in constants_decom...h)
+ ! sets number of nodes per element
ngnod = esize
! reads node coordinates
@@ -388,17 +388,7 @@
close(98)
if( nspec2D_moho > 0 ) print*, ' nspec2D_moho = ', nspec2D_moho
-! Percy & JPA
-
- ! ------------------------------------------------------------
- ! Reading fault elements
- ! -------------------------------------------------------------
call read_fault_files()
-
- !--------------------------------------------------------------
- ! close_fault_crack
- !--------------------------------------------------------------
-
call close_faults(nodes_coords,elmnts,nspec,nnodes,esize)
end subroutine read_mesh_files
@@ -537,22 +527,11 @@
! sup_neighbour, nsize, &
! nparts, part, nfaces_coupled, faces_coupled)
- ! FAULT-Percy , re-partitioning the fault into same partition .
- ! integer :: nfaces_coupled
- ! integer, dimension(:,:), pointer :: faces_coupled
-
- ! re-partioning puts faults coupled elements into same partition.
-! call faultside1_faultside2_repartitioning (nspec, nnodes, elmnts, &
-! count_def_mat, mat(1,:) , &
-! sup_neighbour, nsize, &
-! nparts, part)
-! FAULT : output part(0) : contains all fault elements
-! : fault_ispec1,fault_ispec2 (fault elements side1 and side2)
-! : fault_iface1,fault_iface2 (fault faces side2 and side2)
+ ! move all fault elements to the same partition (proc=0)
call fault_collect_elements(nspec,nnodes,elmnts, &
- sup_neighbour,esize,nsize,nparts,part)
+ sup_neighbour,esize,nsize,nparts,part)
-! re-partitioning puts moho-surface coupled elements into same partition
+ ! re-partitioning puts moho-surface coupled elements into same partition
call moho_surface_repartitioning (nspec, nnodes, elmnts, &
sup_neighbour, nsize, nparts, part, &
nspec2D_moho,ibelm_moho,nodes_ibelm_moho )
@@ -609,20 +588,6 @@
stop 'error file open Database'
endif
-! FAULT : procxxxxx_fault...
-
- ! opens output file
- write(prname, "(i6.6,'_Database_fault')") ipart
- open(unit=16,file=outputpath_name(1:len_trim(outputpath_name))//'/proc'//prname,&
- status='unknown', action='write', form='formatted', iostat = ierr)
- if( ierr /= 0 ) then
- print*,'error file open:',outputpath_name(1:len_trim(outputpath_name))//'/proc'//prname
- print*
- print*,'check if path exists:',outputpath_name(1:len_trim(outputpath_name))
- stop 'error file open Database'
- endif
-
-
! gets number of nodes
call write_glob2loc_nodes_database(15, ipart, nnodes_loc, nodes_coords, &
glob2loc_nodes_nparts, glob2loc_nodes_parts, &
@@ -652,12 +617,6 @@
glob2loc_elmnts, glob2loc_nodes_nparts, &
glob2loc_nodes_parts, glob2loc_nodes, part, mat, ngnod, 2)
-
-! FAULT : Writting out procxxxx_..._fault
- call write_partion_fault_database(16, ipart, nspec, &
- glob2loc_elmnts,part)
-
-
! writes out absorbing/free-surface boundaries
call write_boundaries_database(15, ipart, nspec, nspec2D_xmin, nspec2D_xmax, nspec2D_ymin, &
nspec2D_ymax, nspec2D_bottom, nspec2D_top, &
@@ -689,8 +648,21 @@
nspec2D_moho,ibelm_moho,nodes_ibelm_moho)
close(15)
+
+ ! write fault database
+ write(prname, "(i6.6,'_Database_fault')") ipart
+ open(unit=16,file=outputpath_name(1:len_trim(outputpath_name))//'/proc'//prname,&
+ status='unknown', action='write', form='formatted', iostat = ierr)
+ if( ierr /= 0 ) then
+ print*,'error file open:',outputpath_name(1:len_trim(outputpath_name))//'/proc'//prname
+ print*
+ print*,'check if path exists:',outputpath_name(1:len_trim(outputpath_name))
+ stop
+ endif
+ call write_fault_database(16, ipart, nspec, glob2loc_elmnts,part)
close(16)
+
end do
print*, 'partitions: '
print*, ' num = ',nparts
Modified: seismo/3D/FAULT_SOURCE/branches/new_fault_db/decompose_mesh_SCOTCH/fault_scotch.f90
===================================================================
--- seismo/3D/FAULT_SOURCE/branches/new_fault_db/decompose_mesh_SCOTCH/fault_scotch.f90 2011-05-05 19:52:23 UTC (rev 18322)
+++ seismo/3D/FAULT_SOURCE/branches/new_fault_db/decompose_mesh_SCOTCH/fault_scotch.f90 2011-05-06 00:24:51 UTC (rev 18323)
@@ -5,7 +5,7 @@
type fault_type
private
- integer :: nspec,eta
+ integer :: nspec
integer, dimension(:), pointer :: ispec1, ispec2, iface1, iface2
end type fault_type
@@ -19,21 +19,23 @@
real(kind=CUSTOM_REAL), parameter :: FAULT_GAP_TOLERANCE = 1e-2_CUSTOM_REAL
- public :: read_fault_files, fault_collect_elements, close_faults
+ public :: read_fault_files, fault_collect_elements, close_faults, write_fault_database
CONTAINS
!==========================================================================================
Subroutine read_fault_files
- integer :: nbfaults=0, iflt, ier
+ integer :: nbfaults, iflt, ier
open(101,file='../DATA/FAULT/Par_file_faults.in',status='old',action='read',iostat=ier)
if (ier==0) then
read(101,*) nbfaults
else
+ nbfaults = 0
print *, 'Par_file.in not found: assume no faults'
endif
+ close(101)
if (nbfaults>0) then
allocate(faults(nbfaults))
@@ -42,7 +44,6 @@
enddo
endif
- close(101)
end subroutine read_fault_files
@@ -232,6 +233,8 @@
!jpa: why do this? does it always help balancing ?
!pgb: Yes, bulk elements in processor 0 are taken out and redistributed.
!pgb: leaving more space for fault elements.
+ !jpa: But if the number of fault elements is much smaller than nproc_null
+ ! we will end up with a very UNbalanced proc 0 !
ipart=0
do i = 1, nproc_null
if ( ipart == nproc ) ipart = 0
@@ -239,12 +242,9 @@
part(elem_proc_null(i)) = ipart
end do
-
-! Percy , This is needed to get adjacent element by common face.
call mesh2dual_ncommonnodes(nelmnts, nnodes, nsize, sup_neighbour, &
elmnts, xadj, adjncy, nnodes_elmnts, &
nodes_elmnts, max_neighbour, 4)
-
! coupled elements
! ---------------
@@ -280,7 +280,7 @@
! ---------------------------------------------------------------------------------------------------
! write one block for each fault
- subroutine write_fault_partition_database(IIN_database, iproc, nelmnts, &
+ subroutine write_fault_database(IIN_database, iproc, nelmnts, &
glob2loc_elmnts, part)
! include './constants_decompose_mesh_SCOTCH.h'
@@ -340,7 +340,7 @@
enddo
- end subroutine write_fault_partition_database
+ end subroutine write_fault_database
end module fault_scotch
Modified: seismo/3D/FAULT_SOURCE/branches/new_fault_db/src/create_regions_mesh.f90
===================================================================
--- seismo/3D/FAULT_SOURCE/branches/new_fault_db/src/create_regions_mesh.f90 2011-05-05 19:52:23 UTC (rev 18322)
+++ seismo/3D/FAULT_SOURCE/branches/new_fault_db/src/create_regions_mesh.f90 2011-05-06 00:24:51 UTC (rev 18323)
@@ -53,8 +53,8 @@
! create the different regions of the mesh
use create_regions_mesh_ext_par
- use fault_object, only: fault_read_input, fault_setup, fault_save_arrays_test, fault_save_arrays, &
- fault_db
+ use fault_object, only: fault_read_input, fault_setup, fault_save_arrays, &
+ fault_db !, fault_save_arrays_test
implicit none
@@ -193,7 +193,7 @@
! If the mesh contains faults we split the fault nodes
! and generate the fault database ...
! The node splitting procedure changes ibool size (nglob)
-! and creates Kevin_voigt_eta values .(0 : no damping).
+! and creates Kelvin_voigt_eta values .(0 : no damping).
! crm_ext_setup_indexing : computes xstore , ystore , zstore.
@@ -204,7 +204,6 @@
!NEW : Here loading fault ispec and fault iface.
call fault_read_input()
-
if (allocated(fault_db)) call fault_setup (ibool,xstore,ystore,zstore,nspec,nglob,prname,myrank)
! else
! call crm_ext_setup_indexing(ibool, &
@@ -306,9 +305,7 @@
c55store,c56store,c66store, &
ispec_is_acoustic,ispec_is_elastic,ispec_is_poroelastic)
-!Percy : save fault database
-
- call fault_save_arrays_test(prname,IOUT) ! for debugging
+! call fault_save_arrays_test(prname,IOUT) ! for debugging
call fault_save_arrays(prname,IOUT)
! computes the approximate amount of static memory needed to run the solver
Modified: seismo/3D/FAULT_SOURCE/branches/new_fault_db/src/fault_object.f90
===================================================================
--- seismo/3D/FAULT_SOURCE/branches/new_fault_db/src/fault_object.f90 2011-05-05 19:52:23 UTC (rev 18322)
+++ seismo/3D/FAULT_SOURCE/branches/new_fault_db/src/fault_object.f90 2011-05-06 00:24:51 UTC (rev 18323)
@@ -25,7 +25,7 @@
! This module was written by:
! Percy Galvez, Jean-Paul Ampuero and Tarje Nissen-Meyer
-! Percy : New version with split nodes done in CUBIT.
+! New version with split nodes done in CUBIT.
module fault_object
@@ -76,22 +76,24 @@
!=================================================================================================================
subroutine fault_read_input()
- integer :: nb
+ integer :: nb, i,ier
+ real(kind=CUSTOM_REAL) :: eta
- integer :: i,ier
-
nb = 0
open(unit=100,file='DATA/FAULT/Par_file_faults.in',status='old',action='read',iostat=ier)
if (ier==0) then
- read(100,*) nb
- allocate(fault_db(nb))
- do i=1,nb
- fault_db(i)%eta
- enddo
- else
- write(6,*) 'File Par_file_faults.in does not exist '
- return
+ read(100,*) nb
+ read(100,*) eta
+ if (nb>0) then
+ allocate(fault_db(nb))
+ do i=1,nb
+ fault_db(i)%eta = eta
+ enddo
+ endif
+ else
+ write(6,*) 'File Par_file_faults.in does not exist '
+ return
end if
close(100)
@@ -103,10 +105,6 @@
!==================================================================================================================
subroutine fault_setup(ibool,xstore,ystore,zstore,nspec,nglob,prname,myrank)
-!Percy : mat_ext_mesh(i,ispec) : material index properties
-! Domain tags for each element are in mat_ext_mesh(1,:)
- use generate_databases_par, only : mat_ext_mesh
-
integer, intent(in) :: nspec ! number of spectral elements in each block
double precision, dimension(NGLLX,NGLLY,NGLLZ,nspec), intent(in) :: xstore,ystore,zstore
integer, dimension(NGLLX,NGLLY,NGLLZ,nspec), intent(in) :: ibool
@@ -199,7 +197,6 @@
integer,intent(in) :: myrank
character(len=256), intent(in) :: prname ! 'proc***'
-
integer :: ifault,iflt
integer, dimension(3,NGLLSQUARE,nspec*6) :: ijk1, ijk2
integer :: ijk_face(3,NGLLX,NGLLY)
@@ -215,38 +212,41 @@
do iflt=1,size(faults)
- read(IIN,*) faults(iflt)%nspec
- if (faults(iflt)%nspec == 0) cycle
- allocate(faults(iflt)%ispec1(fdb%nspec))
- allocate(faults(iflt)%iface1(fdb%nspec))
- allocate(faults(iflt)%ijk1(3,NGLLX*NGLLY,fdb%nspec))
+ read(IIN,*) faults(iflt)%nspec
- allocate(faults(iflt)%ispec2(fdb%nspec))
- allocate(faults(iflt)%iface2(fdb%nspec))
- allocate(faults(iflt)%ijk2(3,NGLLX*NGLLY,fdb%nspec))
- ! loading ispec1 ispec2 iface1 iface2 of fault elements.
- do i=1,faults(iflt)%nspec
- read(IIN,*) faults(iflt)%ispec1(i),faults(iflt)%ispec2(i),faults(iflt)%iface1(i),faults(iflt)%iface2(i)
- enddo
+ if (faults(iflt)%nspec == 0) cycle
- do ifault=1,faults(iflt)%nspec
- ! we have identified a new fault element on fault side 1
- iface_ref1 = faults(iflt)%iface1(ifault)
- iface_ref2 = faults(iflt)%iface2(ifault)
- ! gets i,j,k indices of GLL nodes in element face
- call get_element_face_gll_indices(iface_ref1,ijk_face1,NGLLX,NGLLY)
- call get_element_face_gll_indices(iface_ref2,ijk_face2,NGLLX,NGLLY)
- igll = 0
- do j=1,NGLLY
- do i=1,NGLLX
- igll = igll + 1
- ijk1(:,igll,ifault)=ijk_face1(:,i,j) ! saving gll points of side 1 , needed for iulk1.
- ijk2(:,igll,ifault)=ijk_face2(:,i,j) ! saving gll points of side 2 , needed for iulk1.
- enddo
- enddo
+ allocate(faults(iflt)%ispec1(fdb%nspec))
+ allocate(faults(iflt)%iface1(fdb%nspec))
+ allocate(faults(iflt)%ijk1(3,NGLLX*NGLLY,fdb%nspec))
+
+ allocate(faults(iflt)%ispec2(fdb%nspec))
+ allocate(faults(iflt)%iface2(fdb%nspec))
+ allocate(faults(iflt)%ijk2(3,NGLLX*NGLLY,fdb%nspec))
+
+ ! loading ispec1 ispec2 iface1 iface2 of fault elements.
+ do i=1,faults(iflt)%nspec
+ read(IIN,*) faults(iflt)%ispec1(i),faults(iflt)%ispec2(i),faults(iflt)%iface1(i),faults(iflt)%iface2(i)
+ enddo
+
+ do ifault=1,faults(iflt)%nspec
+ ! we have identified a new fault element on fault side 1
+ iface_ref1 = faults(iflt)%iface1(ifault)
+ iface_ref2 = faults(iflt)%iface2(ifault)
+ ! gets i,j,k indices of GLL nodes in element face
+ call get_element_face_gll_indices(iface_ref1,ijk_face1,NGLLX,NGLLY)
+ call get_element_face_gll_indices(iface_ref2,ijk_face2,NGLLX,NGLLY)
+ igll = 0
+ do j=1,NGLLY
+ do i=1,NGLLX
+ igll = igll + 1
+ ijk1(:,igll,ifault)=ijk_face1(:,i,j) ! saving gll points of side 1 , needed for ibulk1
+ ijk2(:,igll,ifault)=ijk_face2(:,i,j) ! saving gll points of side 2 , needed for ibulk2
+ enddo
enddo
- faults(iflt)%ijk1 = ijk1(:,:,1:faults(iflt)%nspec)
- faults(iflt)%ijk2 = ijk2(:,:,1:faults(iflt)%nspec)
+ enddo
+ faults(iflt)%ijk1 = ijk1(:,:,1:faults(iflt)%nspec)
+ faults(iflt)%ijk2 = ijk2(:,:,1:faults(iflt)%nspec)
enddo
More information about the CIG-COMMITS
mailing list