[cig-commits] r20825 - seismo/3D/FAULT_SOURCE/branches/new_fault_db/src

[surendra at geodynamics.org]

Author: surendra
Date: 2012-10-11 22:44:56 -0700 (Thu, 11 Oct 2012)
New Revision: 20825

Modified:
   seismo/3D/FAULT_SOURCE/branches/new_fault_db/src/fault_solver.f90
Log:
Benchmarked rate-and-state for ageing law TPV102.  Results now match with those of other modelers for TPV102

Modified: seismo/3D/FAULT_SOURCE/branches/new_fault_db/src/fault_solver.f90
===================================================================
--- seismo/3D/FAULT_SOURCE/branches/new_fault_db/src/fault_solver.f90	2012-10-11 19:52:31 UTC (rev 20824)
+++ seismo/3D/FAULT_SOURCE/branches/new_fault_db/src/fault_solver.f90	2012-10-12 05:44:56 UTC (rev 20825)
@@ -558,16 +558,16 @@
       call init_2d_distribution(bc%rsf%C  ,bc%coord,IIN_PAR,nC)
       call init_2d_distribution(bc%rsf%T  ,bc%coord,IIN_PAR,nForcedRup)
 
-    ! WARNING : Not general enough
-    nglob_bulk = size(vel,2)
-    allocate(init_vel(3,nglob_bulk))
-    init_vel = 0._CUSTOM_REAL
-    init_vel(1,bc%ibulk1) =  -bc%rsf%V_init/2._CUSTOM_REAL
-    init_vel(1,bc%ibulk2) =  bc%rsf%V_init/2._CUSTOM_REAL
-    where(ystore > 0) init_vel(1,:) = -V_init/2._CUSTOM_REAL
-    where(ystore < 0) init_vel(1,:) = V_init/2._CUSTOM_REAL
-    !init_vel = rotate(bc,init_vel,-1) ! directly assigned in global coordinates here as it is a simplified case
-    vel = vel + init_vel
+!!$    ! WARNING : Not general enough
+!!$    nglob_bulk = size(vel,2)
+!!$    allocate(init_vel(3,nglob_bulk))
+!!$    init_vel = 0._CUSTOM_REAL
+!!$    init_vel(1,bc%ibulk1) =  -bc%rsf%V_init/2._CUSTOM_REAL
+!!$    init_vel(1,bc%ibulk2) =  bc%rsf%V_init/2._CUSTOM_REAL
+!!$    where(ystore > 0) init_vel(1,:) = -V_init/2._CUSTOM_REAL
+!!$    where(ystore < 0) init_vel(1,:) = V_init/2._CUSTOM_REAL
+!!$    !init_vel = rotate(bc,init_vel,-1) ! directly assigned in global coordinates here as it is a simplified case
+!!$    vel = vel + init_vel
     
       W1=15000._CUSTOM_REAL
       W2=7500._CUSTOM_REAL
@@ -922,7 +922,8 @@
     
     ! Subtract initial stress
     T = T - bc%T0
-    
+    if(RATE_AND_STATE) T(1,:) = T(1,:) + TxExt
+
     ! Update slip acceleration da=da_free-T/(0.5*dt*Z)
     dA(1,:) = dA(1,:) - T(1,:)/(bc%Z*half_dt)
     dA(2,:) = dA(2,:) - T(2,:)/(bc%Z*half_dt)