fixed bug in hh_visco.c

fixed bug which caused segmentation fault by using WRITE_MODELFILE=2

src/hh_visco.c

@@ -54,11 +54,12 @@ void model_visco(float  ***  rho, float ***  pi, float ***  u,
 	/*-----------------material property definition -------------------------*/
 
 	/* parameters for layer 1 */
-	const float vp1=3500.0, vs1=2000.0, rho1=2000.0, h=1000000, qp1=20.0, qs1=10.0;
+	const float vp1=3500.0, vs1=2000.0, rho1=2000.0, h=100000.0, qp1=20.0, qs1=10.0;
 
 	/* parameters for layer 2 */
 	const float vp2=5700.0, vs2=3400.0, rho2=2500.0, qp2=200.0, qs2=100.0;
-
+	
+	
 
 	/*internal switch for writing all models to file (WRITE_MODELFILES=1)
 	 * or just density (WRITE_MODELFILES=0)
@@ -88,12 +89,12 @@ void model_visco(float  ***  rho, float ***  pi, float ***  u,
 		/* in the viscoelastic case : reference frequency where no velocity dispersion occurs.
 		 * if FREF is not given in input file, the largest center source frequency FC
 		 * as specified in input file is used (not that the relation : FC=1/TS) is used here)*/
-
+		
 		if (FREF==0.0) ws=2.0*PI/TS;
 		else ws=2.0*PI*FREF;
 
 		//fprintf(FP,"MYID=%d \t\t ws=%5.5f \t pts=%5.5f \t FL=%5.5f \n ",MYID,ws,pts[l],FL[l]);
-
+		
 		/* loop over global grid */
 		for (k=1;k<=NZG;k++){
 			for (i=1;i<=NXG;i++){
@@ -102,7 +103,7 @@ void model_visco(float  ***  rho, float ***  pi, float ***  u,
 					 * modify below this point for VISCOELASTIC model definition
 					 *=========================================================
 					 */
-
+					
 					/*note that "y" is used for the vertical coordinate*/
 					/* calculate vertical coordinate in m */
 
@@ -155,7 +156,7 @@ void model_visco(float  ***  rho, float ***  pi, float ***  u,
 							taup[jj][ii][kk]=TAU;
 						}
 
-						if (WRITE_MODELFILES) {
+						if (WRITE_MODELFILES==1) {
 							pwavemod[jj][ii][kk]=Vp;
 							swavemod[jj][ii][kk]=Vs;
 							qsmod[jj][ii][kk]=2/taus[jj][ii][kk];
@@ -163,13 +164,14 @@ void model_visco(float  ***  rho, float ***  pi, float ***  u,
 						}
 
 					}
-				}
+				} 
 			}
 		}
+		
 	}
 
 	/* each PE writes his model to disk */
-
+	
 	/* all models are written to file */
 	if (WRITE_MODELFILES==1) {
 		sprintf(modfile,"%s.SOFI3D.pi",MFILE);
@@ -219,7 +221,7 @@ void model_visco(float  ***  rho, float ***  pi, float ***  u,
 	}
 
 	free_vector(pts,1,L);
-	if (WRITE_MODELFILES) {
+	if (WRITE_MODELFILES==1) {
 		free_f3tensor(pwavemod,0,NY+1,0,NX+1,0,NZ+1);
 		free_f3tensor(swavemod,0,NY+1,0,NX+1,0,NZ+1);
 		if ((L) && (TAU==0.0)){
@@ -227,5 +229,5 @@ void model_visco(float  ***  rho, float ***  pi, float ***  u,
 			free_f3tensor(qsmod,0,NY+1,0,NX+1,0,NZ+1);
 		}
 	}
-
+	
 }