Ticket #38522: relax.diff

src/export.f90

@@ -920 +920 @@
        WRITE (15,'("#        x1         x2         x3          yr        yz", &
                  & "       slip strike-slip  dip-slip   velocity     ss vel", &
                  & "     ds vel       taus      sig11      sig12      sig13      sig22      sig23      sig33")')
-
-       WRITE (15,'(18ES11.3E2)') ((n(k)%patch(i1,i2)%x1,n(k)%patch(i1,i2)%x2,n(k)%patch(i1,i2)%x3, &
+       WRITE (15,'(18ES11.3E2)') ((n(k)%patch(i1,i2)%x1,n(k)%patch(i1,i2)%x3,n(k)%patch(i1,i2)%x3, &
                                   n(k)%patch(i1,i2)%lx,n(k)%patch(i1,i2)%lz, &
                                   n(k)%patch(i1,i2)%slip, &
                                   n(k)%patch(i1,i2)%ss, &

src/friction3d.f90

@@ -312 +312 @@
           patch(j2,j3)%x3=x3
           patch(j2,j3)%lx=yr
           patch(j2,j3)%lz=zr
+          patch(j2,j3)%sig=s
 
           ! discard out-of-bound locations
           IF (  (x1 .GT. DBLE(sx1/2-1)*dx1) .OR. (x1 .LT. -DBLE(sx1/2)*dx1) &
@@ -329 +330 @@
        
           ! traction = sigma . n
           s=sig(i1,i2,i3)
-          patch(j2,j3)%sig=s
           t=s .tdot. n
 
           ! signed normal component
@@ -426 +426 @@
        DO j3=1,px3
           ! loop in the strike direction
           DO j2=1,px2
+             ! cumulative creep
+             n(k)%patch(j2,j3)%slip=n(k)%patch(j2,j3)%slip+dt*n(k)%patch(j2,j3)%v
+
              ! cumulative strike-direction creep
              n(k)%patch(j2,j3)%ss=n(k)%patch(j2,j3)%ss+dt*n(k)%patch(j2,j3)%vss
 
              ! cumulative dip-direction creep
              n(k)%patch(j2,j3)%ds=n(k)%patch(j2,j3)%ds+dt*n(k)%patch(j2,j3)%vds
-
-             ! cumulative creep
-             n(k)%patch(j2,j3)%slip=(n(k)%patch(j2,j3)%ds**2+n(k)%patch(j2,j3)%ss**2)**0.5
           END DO
        END DO
 

src/input.f90

@@ -818 +818 @@
 
              ALLOCATE(in%n(k)%patch(in%n(k)%px2,in%n(k)%px3),STAT=iostatus)
              IF (iostatus>0) STOP "could not allocate the fault patches"
-             in%n(k)%patch(:,:)=SLIPPATCH_STRUCT(0,0,0,0,0,0,0,0,0,0,0,0,TENSOR(0,0,0,0,0,0))
 
 #ifdef VTK
              ! export the afterslip segment in VTK format
@@ -860 +859 @@
                in%inter%s(k)%x,in%inter%s(k)%y,in%inter%s(k)%z, &
                in%inter%s(k)%length,in%inter%s(k)%width, &
                in%inter%s(k)%strike,in%inter%s(k)%dip,in%inter%s(k)%rake
-          in%inter%s(k)%opening=0
 
           ! copy the input format for display
           in%inter%sc(k)=in%inter%s(k)
@@ -1025 +1023 @@
                   in%events(e)%s(k)%length,in%events(e)%s(k)%width, &
                   in%events(e)%s(k)%strike,in%events(e)%s(k)%dip,in%events(e)%s(k)%rake, &
                   in%events(e)%s(k)%beta
-             in%events(e)%s(k)%opening=0
 
              SELECT CASE(iostatus)
              CASE (1:)

src/relax.f90

@@ -270 +270 @@
   IF (in%isoutputrelax) THEN
      ALLOCATE(inter1(in%sx1+2,in%sx2,2),inter2(in%sx1+2,in%sx2,2),inter3(in%sx1+2,in%sx2,2),STAT=iostatus)
      IF (iostatus>0) STOP "could not allocate memory for postseismic displacement"
-     inter1=0;inter2=0;inter3=0;
   END IF
 
   v1=0;v2=0;v3=0;u1=0;u2=0;u3=0;gamma=0;t1=0;t2=0;t3=0
@@ -333 +332 @@
                          in%dx1,in%dx2,in%dx3,in%sx1,in%sx2,in%sx3/2,in%wdir,0)
   
   ! export equivalent body forces
-  IF (isoutput(in%skip,t,0,in%odt,oi,in%events(e)%time)) THEN
+  IF (isoutput(in%skip,t,i,in%odt,oi,in%events(e)%time)) THEN
 #ifdef GRD_EQBF
      IF (in%isoutputgrd) THEN
         CALL exportgrd(v1,v2,v3,in%sx1,in%sx2,in%sx3/2, &

wscript

@@ -53 +53 @@
 def configure(cnf):
     cnf.load('compiler_c compiler_fc')
     
-    # # Check for -zero or -finit-local-zero.  We check it first, before
-    # # setting CFLAGS, because optimization can ruin this test.
-    # frag="program main\n  INTEGER :: foo\n  call exit(foo)\n" + "end program main\n"
-    # zero_flags=['-zero','-finit-local-zero']
-    # if cnf.options.zero_flag:
-    #     zero_flags=[cnf.options.zero_flag]
-    # found_zero=False
-    # for flag in zero_flags:
-    #     try:
-    #         cnf.check_fc(fragment=frag,msg="Checking zero option " + flag,
-    #                      fcflags=flag,uselib_store='zero',execute=True)
-    #     except:
-    #         continue
-    #     else:
-    #         found_zero=True
-    #         break
-    # if not found_zero:
-    #     cnf.fatal("Could not find an option for zero'ing initial values")
+    # Check for -zero or -finit-local-zero.  We check it first, before
+    # setting CFLAGS, because optimization can ruin this test.
+    frag="program main\n  INTEGER :: foo\n  call exit(foo)\n" + "end program main\n"
+    zero_flags=['-zero','-finit-local-zero']
+    if cnf.options.zero_flag:
+        zero_flags=[cnf.options.zero_flag]
+    found_zero=False
+    for flag in zero_flags:
+        try:
+            cnf.check_fc(fragment=frag,msg="Checking zero option " + flag,
+                         fcflags=flag,uselib_store='zero',execute=True)
+        except:
+            continue
+        else:
+            found_zero=True
+            break
+    if not found_zero:
+        cnf.fatal("Could not find an option for zero'ing initial values")
 
     # We set the flags here 
     if not cnf.env.CFLAGS:
@@ -234 +234 @@
                         'src/input.f90',
                         'src/mkl_dfti.f90'],
                 includes=['build'],
-                use=['gmt','proj','openmp','fftw','imkl','cpp','length'],
+                use=['gmt','proj','openmp','fftw','imkl','zero','cpp','length'],
                 target='relax'
                 )