Commit ba2301c1 authored by tilman.metz's avatar tilman.metz

Small Bugfix

At some computers rounding errors occur
The coefficient alpha_prime should be > 0
But somehow the calculation can end up negativ for the coefficient near to the pml
border (should be zero). If that happens it will be set to zero now.
parent 0a9815f5
......@@ -133,7 +133,10 @@ void CPML_coeff(st_pml_coeff *pml_coeff)
pml_coeff->K_x[i] = 1.0 + (k_max_PML - 1.0) * pow(position_norm,npower);
pml_coeff->alpha_prime_x[i] = alpha_max_PML * (1.0 - position_norm);
if(pml_coeff->alpha_prime_x[i] < 0.0){ fprintf(FP,"ERROR:alpha_prime_x[i] < 0.0, i %d", i);}
if(pml_coeff->alpha_prime_x[i] < 0.0)
{ fprintf(FP,"ERROR:alpha_prime_x[i] < 0.0, alpha_prime=%e, i %d", pml_coeff->alpha_prime_x[i],i);
if (i==1) pml_coeff->alpha_prime_x[i]=0.0;
}
pml_coeff->b_x[i] = exp(- (d_x[i] / pml_coeff->K_x[i] + pml_coeff->alpha_prime_x[i]) * DT);
......@@ -202,8 +205,11 @@ void CPML_coeff(st_pml_coeff *pml_coeff)
pml_coeff->alpha_prime_y[i] = alpha_max_PML * (1.0 - position_norm);
/* just in case, for -5 at the end */
if(pml_coeff->alpha_prime_y[i] < 0.0){ fprintf(FP,"ERROR:alpha_prime_y[i] < 0.0, i %d", i);}
if(pml_coeff->alpha_prime_y[i] < 0.0){ fprintf(FP,"ERROR:alpha_prime_y[i] < 0.0,alpha_prime_y[i] %e, i %d\n", pml_coeff->alpha_prime_y[1],i);
if (i==1) pml_coeff->alpha_prime_y[i]=0.0;
}
pml_coeff->b_y[i] = exp(- (d_y[i] / pml_coeff->K_y[i] + pml_coeff->alpha_prime_y[i]) * DT);
/* avoid division by zero outside the PML */
......@@ -266,7 +272,9 @@ void CPML_coeff(st_pml_coeff *pml_coeff)
pml_coeff->alpha_prime_z[i] = alpha_max_PML * (1.0 - position_norm);
/* just in case, for -5 at the end */
if(pml_coeff->alpha_prime_z[i] < 0.0){ fprintf(FP,"ERROR:alpha_prime_z[i] < 0.0, i %d", i);}
if(pml_coeff->alpha_prime_z[i] < 0.0){ fprintf(FP,"ERROR:alpha_prime_z[i] < 0.0,alpha_prime_z[i]==%e, i %d\n",pml_coeff->alpha_prime_z[i], i);
if (i==1) pml_coeff->alpha_prime_z[i]=0.0;
}
pml_coeff->b_z[i] = exp(- (d_z[i] / pml_coeff->K_z[i] + pml_coeff->alpha_prime_z[i]) * DT);
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment