c*********************************************************************** c---+----1----+----2----+----3----+----4----+----5----+----6----+----7-- c*********************************************************************** !--------------------------------------------------------------------------------------------------- ! BOP ! ! ROUTINE: DERIVS ! ! INTERFACE: SUBROUTINE DERIVS ( & t, & y, ! supplied - array of stresses and velocities on each cell & dydt) ! returned - array of derivatives of stresses and velocities ! ! USES: USE Green_sshr ! Make data in module "Green_sshr" available use cblock ! use data and variables in cblock ! ! RETURN VALUE: c implicit none ! Not yet converted to using "implicit none" and explicitly typing all variables c :: ! Description implicit real*8 (a-h,o-z) common /newderiv/ p1(4), nnn, myid, umm(4) dimension y(1),dydt(1) ! ! NECESSARY PRE-CONDITIONS: ! none ! ! SIDE EFFECTS: ! none ! ! POST_CONDITIONS UPON EXIT ! none ! ! DESCRIPTION: ! Called by MAIN, by RKQC and by RK4 ! ! ! ! calc derivatives (rt sides) of ode's t isn't in ode's ! shear stress - sshr(src,obs) ! dimension y,dydt by 2*nf in main; eqns for non-bc cells only ! vpl,vel,u,tauss always >0 ! calc dtau/dt from disloc stresses, then dv/dt from fault law ! ! uses radiation damping (quasi-dynamic) approximation to full elastodynamics ! uses slip law version of evolution equation for rate and state friction ! ! EOP !--------------------------------------------------------------------------------------------------- e10=p1(1) xldis=p1(2) vpl=p1(3) xmu=p1(4) ne=nnn c ne=2 c print *, myid, nnn c print 5, (sxc(i), szc(i), skey1(i), myid, i=1,ne) 5 format (2(e10.4, 1x) , i12, 1x, i6) rigid = 3.0d4 ! rigidity in MPa --- NK vs = 3.0d3 ! shear-wave speed in m/s --- NK rigid = rigid * 10.0d0 ! MPa ---> bar NK vs = vs * 1.0d3 * 60.0d0 * 60.0d0 * 24.0d0 * 365.25d0 ! mm/year NK c print *, myid, ' calling buildtree ne = ', ne CALL buildtree(ne, sxc(1), szc(1), sdx(1), sdz(1), + y(1+ne), skey1(1), umm(1)) do 10 n=1,ne ! obs pt: y(n) is tau, y(n+ne) is vel strs = stauvp(n) ! only bc cells slipping k = n + ne ! radiation damping denom = 2.0d0 * vs * sca( n ) + rigid * y( k ) ! radiation damping coefv = 2.0d0 * vs * y( k ) /denom ! radiation damping coeft1 = 2.0d0 * vs * sca( n ) / denom ! radiation damping coeft2 = rigid * y( k ) / denom ! radiation damping CALL sumtree(sxc(n), szc(n), strs) strs = 2.*xmu*strs ! this done in GREENS_FUNCTION for noMP version; not done in tree.c arg = vpl/y(k)+e10 stauss(n) = -scab(n)* log(arg) trm = (y(n)-stauss(n))* y(k)/xldis dydt( n ) = coeft1 * strs - coeft2 * trm ! dtau/dt radiation damping dydt( k ) = coefv * ( strs + trm ) ! dvel/dt radiation damping c print *, myid, y(n), y(k), dydt(n), dydt(k), skey1(n) 10 continue CALL forgettree ! Multipole version c print *, 'tree forgot processor ', myid return end c*********************************************************************** c---+----1----+----2----+----3----+----4----+----5----+----6----+----7-- c***********************************************************************