[cig-commits] r17878 - seismo/2D/SPECFEM2D/trunk

[dkomati1 at geodynamics.org]

Author: dkomati1
Date: 2011-02-15 16:29:33 -0800 (Tue, 15 Feb 2011)
New Revision: 17878

Modified:
   seismo/2D/SPECFEM2D/trunk/part_unstruct.F90
Log:
fixed a bug in subroutine read_abs_surface, which reads an absorbing surface: the number of nodes that form the absorbing edge currently must always be equal to two, even in the case of 9-node elements


Modified: seismo/2D/SPECFEM2D/trunk/part_unstruct.F90
===================================================================
--- seismo/2D/SPECFEM2D/trunk/part_unstruct.F90	2011-02-15 23:45:21 UTC (rev 17877)
+++ seismo/2D/SPECFEM2D/trunk/part_unstruct.F90	2011-02-16 00:29:33 UTC (rev 17878)
@@ -275,8 +275,9 @@
 
   !-----------------------------------------------
   ! Read absorbing surface.
-  ! 'abs_surface' contains 1/ element number, 2/ number of nodes that form the abs surface,
-  ! 3/ first node on the abs surface, 4/ second node on the abs surface, if relevant (if 2/ is equal to 2)
+  ! 'abs_surface' contains 1/ element number, 2/ number of nodes that form the absorbing edge
+  ! (which currently must always be equal to two, see comment below),
+  ! 3/ first node on the abs surface, 4/ second node on the abs surface
   !-----------------------------------------------
   subroutine read_abs_surface(filename, num_start)
 
@@ -286,7 +287,6 @@
   character(len=256), intent(in)  :: filename
   integer, intent(in)  :: num_start
 
-
   integer  :: i,ier
 
   open(unit=994, file=trim(filename), form='formatted' , status='old', action='read', iostat=ier)
@@ -301,7 +301,14 @@
 
   do i = 1, nelemabs
      read(994,*) abs_surface(1,i), abs_surface(2,i), abs_surface(3,i), abs_surface(4,i)
-
+! The input format is currently limited: only two nodes per element can be listed.
+! If one of your elements has more than one edge along a given absorbing contour (e.g., if that contour has a corner)
+! then list it twice, putting the first edge on the first line and the second edge on the second line.
+! if one of your elements has a single point along the absording contour rather than a full edge, do NOT list it
+! (it would have no weight in the contour integral anyway because it would consist of a single point).
+! If you are using 9-node elements, list only the first and last point of the edge and not the intermediate point
+! located around the middle of the edge; the right 9-node curvature will be restored automatically by the code.
+     if (abs_surface(2,i) /= 2) stop 'only two nodes per element should be listed for absorbing edges'
   enddo
 
   close(994)