Skip to content

Commit 91264fd

Browse files
authored
Merge pull request #6072 from gassmoeller/qualify_sqrt
Qualify std functions everywhere
2 parents 6e2d242 + c002d80 commit 91264fd

File tree

54 files changed

+269
-269
lines changed

Some content is hidden

Large Commits have some content hidden by default. Use the searchbox below for content that may be hidden.

54 files changed

+269
-269
lines changed

benchmarks/advection_in_annulus/advection_in_annulus.cc

Lines changed: 5 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -64,10 +64,10 @@ namespace aspect
6464
const double theta = std::atan2(y,x);
6565
const double f_r = A*r + B/r;
6666
const double g_r = A*r/2 + B*std::log(r)/r + C/r;
67-
const double v_r = g_r*k*sin(k*theta);
68-
const double v_theta = f_r*cos(k*theta);
69-
const double v_x = cos(theta)*v_r - sin(theta)*v_theta;
70-
const double v_y = sin(theta)*v_r + cos(theta)*v_theta;
67+
const double v_r = g_r*k*std::sin(k*theta);
68+
const double v_theta = f_r*std::cos(k*theta);
69+
const double v_x = std::cos(theta)*v_r - std::sin(theta)*v_theta;
70+
const double v_y = std::sin(theta)*v_r + std::cos(theta)*v_theta;
7171
return Point<2> (v_x,v_y);
7272
}
7373

@@ -82,7 +82,7 @@ namespace aspect
8282
const double f_r = 2*r + B/r;
8383
const double g_r = A*r/2 + B*std::log(r)/r + C/r;
8484
const double h_r=(2*g_r-f_r)/r;
85-
return k*h_r*sin(k*theta)+rho_0*gravity*(outer_radius-r);
85+
return k*h_r*std::sin(k*theta)+rho_0*gravity*(outer_radius-r);
8686
}
8787

8888
template <int dim>

benchmarks/annulus/plugin/annulus.cc

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -219,7 +219,7 @@ namespace aspect
219219
const double r = spherical_position[0];
220220
const double theta = spherical_position[1];
221221

222-
const double forcing_term = m(r,k)*k*sin(k*(theta-phase(t))) + rho_0;
222+
const double forcing_term = m(r,k)*k*std::sin(k*(theta-phase(t))) + rho_0;
223223
gravity_magnitude = forcing_term / Annulus_density(pos,k,t, transient);
224224
}
225225

benchmarks/burstedde/burstedde.cc

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -192,7 +192,7 @@ namespace aspect
192192
const double x=p[0];
193193
const double y=p[1];
194194
const double z=p[2];
195-
const double mu=exp(1. - beta * (x*(1.-x)+y*(1.-y) + z*(1.-z)));
195+
const double mu=std::exp(1. - beta * (x*(1.-x)+y*(1.-y) + z*(1.-z)));
196196

197197
out.viscosities[i] = mu;
198198
out.thermal_conductivities[i] = 0.0;
@@ -352,7 +352,7 @@ namespace aspect
352352
const double x=pos[0];
353353
const double y=pos[1];
354354
const double z=pos[2];
355-
const double mu=exp(1. - beta * (x*(1.-x)+y*(1.-y) + z*(1.-z)));
355+
const double mu=std::exp(1. - beta * (x*(1.-x)+y*(1.-y) + z*(1.-z)));
356356

357357
const double dmudx=-beta*(1.-2.*x)*mu;
358358
const double dmudy=-beta*(1.-2.*y)*mu;

benchmarks/crameri_et_al/case_1/crameri_benchmark_1.cc

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -46,7 +46,7 @@ namespace aspect
4646
/**
4747
* Generate a coarse mesh for the geometry described by this class.
4848
* Makes perturbs the top boundary of the box with a function
49-
* of the form z' = amplitude * cos(order * x )
49+
* of the form z' = amplitude * std::cos(order * x )
5050
*/
5151
void create_coarse_mesh (parallel::distributed::Triangulation<dim> &coarse_grid) const override;
5252

benchmarks/hollow_sphere/hollow_sphere.cc

Lines changed: 36 additions & 36 deletions
Original file line numberDiff line numberDiff line change
@@ -73,25 +73,25 @@ namespace aspect
7373

7474
if (mmm == -1)
7575
{
76-
alpha=-gammma*(pow(R2,3)-pow(R1,3))/(pow(R2,3)*log(R1)-pow(R1,3)*log(R2));
77-
beta=-3*gammma*(log(R2)-log(R1))/(pow(R1,3)*log(R2)-pow(R2,3)*log(R1)) ;
76+
alpha=-gammma*(std::pow(R2,3)-std::pow(R1,3))/(std::pow(R2,3)*std::log(R1)-std::pow(R1,3)*std::log(R2));
77+
beta=-3*gammma*(std::log(R2)-std::log(R1))/(std::pow(R1,3)*std::log(R2)-std::pow(R2,3)*std::log(R1)) ;
7878
fr=alpha/(r*r)+beta*r;
79-
gr=-2/(r*r)*(alpha*log(r)+beta/3*pow(r,3)+gammma);
79+
gr=-2/(r*r)*(alpha*std::log(r)+beta/3*std::pow(r,3)+gammma);
8080
}
8181
else
8282
{
83-
alpha=gammma*(mmm+1)*(pow(R1,-3)-pow(R2,-3))/(pow(R1,-mmm-4)-pow(R2,-mmm-4));
84-
beta=-3*gammma*(pow(R1,mmm+1)-pow(R2,mmm+1))/(pow(R1,mmm+4)-pow(R2,mmm+4));
85-
fr=alpha/pow(r,mmm+3)+beta*r;
86-
gr=-2/(r*r)*(-alpha/(mmm+1)*pow(r,-mmm-1)+beta/3*pow(r,3)+gammma);
83+
alpha=gammma*(mmm+1)*(std::pow(R1,-3)-std::pow(R2,-3))/(std::pow(R1,-mmm-4)-std::pow(R2,-mmm-4));
84+
beta=-3*gammma*(std::pow(R1,mmm+1)-std::pow(R2,mmm+1))/(std::pow(R1,mmm+4)-std::pow(R2,mmm+4));
85+
fr=alpha/std::pow(r,mmm+3)+beta*r;
86+
gr=-2/(r*r)*(-alpha/(mmm+1)*std::pow(r,-mmm-1)+beta/3*std::pow(r,3)+gammma);
8787
}
8888

89-
const double v_r =gr*cos(theta);
90-
const double v_theta=fr*sin(theta);
91-
const double v_phi =fr*sin(theta);
92-
const double v_x=sin(theta)*cos(phi)*v_r + cos(theta)*cos(phi)*v_theta-sin(phi)*v_phi;
93-
const double v_y=sin(theta)*sin(phi)*v_r + cos(theta)*sin(phi)*v_theta+cos(phi)*v_phi;
94-
const double v_z=cos(theta)*v_r - sin(theta)*v_theta;
89+
const double v_r =gr*std::cos(theta);
90+
const double v_theta=fr*std::sin(theta);
91+
const double v_phi =fr*std::sin(theta);
92+
const double v_x=std::sin(theta)*std::cos(phi)*v_r + std::cos(theta)*std::cos(phi)*v_theta-std::sin(phi)*v_phi;
93+
const double v_y=std::sin(theta)*std::sin(phi)*v_r + std::cos(theta)*std::sin(phi)*v_theta+std::cos(phi)*v_phi;
94+
const double v_z=std::cos(theta)*v_r - std::sin(theta)*v_theta;
9595

9696
// create a Point<3> (because it has a constructor that takes
9797
// three doubles) and return it (it automatically converts to
@@ -115,21 +115,21 @@ namespace aspect
115115
if (mmm == -1)
116116
{
117117
mur=mu0;
118-
alpha=-gammma*(pow(R2,3)-pow(R1,3))/(pow(R2,3)*log(R1)-pow(R1,3)*log(R2));
119-
beta=-3*gammma*(log(R2)-log(R1))/(pow(R1,3)*log(R2)-pow(R2,3)*log(R1)) ;
120-
gr=-2/(r*r)*(alpha*log(r)+beta/3*pow(r,3)+gammma);
118+
alpha=-gammma*(std::pow(R2,3)-std::pow(R1,3))/(std::pow(R2,3)*std::log(R1)-std::pow(R1,3)*std::log(R2));
119+
beta=-3*gammma*(std::log(R2)-std::log(R1))/(std::pow(R1,3)*std::log(R2)-std::pow(R2,3)*std::log(R1)) ;
120+
gr=-2/(r*r)*(alpha*std::log(r)+beta/3*std::pow(r,3)+gammma);
121121
hr=2/r*gr*mur;
122122
}
123123
else
124124
{
125-
mur=mu0*pow(r,mmm+1);
126-
alpha=gammma*(mmm+1)*(pow(R1,-3)-pow(R2,-3))/(pow(R1,-mmm-4)-pow(R2,-mmm-4));
127-
beta=-3*gammma*(pow(R1,mmm+1)-pow(R2,mmm+1))/(pow(R1,mmm+4)-pow(R2,mmm+4));
128-
gr=-2/(r*r)*(-alpha/(mmm+1)*pow(r,-mmm-1)+beta/3*pow(r,3)+gammma);
125+
mur=mu0*std::pow(r,mmm+1);
126+
alpha=gammma*(mmm+1)*(std::pow(R1,-3)-std::pow(R2,-3))/(std::pow(R1,-mmm-4)-std::pow(R2,-mmm-4));
127+
beta=-3*gammma*(std::pow(R1,mmm+1)-std::pow(R2,mmm+1))/(std::pow(R1,mmm+4)-std::pow(R2,mmm+4));
128+
gr=-2/(r*r)*(-alpha/(mmm+1)*std::pow(r,-mmm-1)+beta/3*std::pow(r,3)+gammma);
129129
hr=(mmm+3)/r*gr*mur;
130130
}
131131

132-
return hr*cos(theta) + rho_0 * gravity * (R2 - r);
132+
return hr*std::cos(theta) + rho_0 * gravity * (R2 - r);
133133
}
134134

135135
template <int dim>
@@ -147,20 +147,20 @@ namespace aspect
147147

148148
if (mmm == -1)
149149
{
150-
alpha=-gammma*(pow(R2,3)-pow(R1,3))/(pow(R2,3)*log(R1)-pow(R1,3)*log(R2));
151-
beta=-3*gammma*(log(R2)-log(R1))/(pow(R1,3)*log(R2)-pow(R2,3)*log(R1)) ;
150+
alpha=-gammma*(std::pow(R2,3)-std::pow(R1,3))/(std::pow(R2,3)*std::log(R1)-std::pow(R1,3)*std::log(R2));
151+
beta=-3*gammma*(std::log(R2)-std::log(R1))/(std::pow(R1,3)*std::log(R2)-std::pow(R2,3)*std::log(R1)) ;
152152
fr=alpha/(r*r)+beta*r;
153-
gr=-2/(r*r)*(alpha*log(r)+beta/3*pow(r,3)+gammma);
153+
gr=-2/(r*r)*(alpha*std::log(r)+beta/3*std::pow(r,3)+gammma);
154154
}
155155
else
156156
{
157-
alpha=gammma*(mmm+1)*(pow(R1,-3)-pow(R2,-3))/(pow(R1,-mmm-4)-pow(R2,-mmm-4));
158-
beta=-3*gammma*(pow(R1,mmm+1)-pow(R2,mmm+1))/(pow(R1,mmm+4)-pow(R2,mmm+4));
159-
fr=alpha/pow(r,mmm+3)+beta*r;
160-
gr=-2/(r*r)*(-alpha/(mmm+1)*pow(r,-mmm-1)+beta/3*pow(r,3)+gammma);
157+
alpha=gammma*(mmm+1)*(std::pow(R1,-3)-std::pow(R2,-3))/(std::pow(R1,-mmm-4)-std::pow(R2,-mmm-4));
158+
beta=-3*gammma*(std::pow(R1,mmm+1)-std::pow(R2,mmm+1))/(std::pow(R1,mmm+4)-std::pow(R2,mmm+4));
159+
fr=alpha/std::pow(r,mmm+3)+beta*r;
160+
gr=-2/(r*r)*(-alpha/(mmm+1)*std::pow(r,-mmm-1)+beta/3*std::pow(r,3)+gammma);
161161
}
162162

163-
return -(6.*gr + 4.*fr) * cos(theta) * mu0 / r;
163+
return -(6.*gr + 4.*fr) * std::cos(theta) * mu0 / r;
164164
}
165165

166166

@@ -366,7 +366,7 @@ namespace aspect
366366
const Point<dim> &pos = in.position[i];
367367
const std::array<double,dim> spos = aspect::Utilities::Coordinates::cartesian_to_spherical_coordinates(pos);
368368
const double r = spos[0];
369-
const double mu = pow(r,mmm+1);
369+
const double mu = std::pow(r,mmm+1);
370370
out.viscosities[i] = mu;
371371

372372
const double theta=spos[2];
@@ -380,15 +380,15 @@ namespace aspect
380380

381381
if (mmm == -1)
382382
{
383-
alpha = -gammma*(pow(R2,3)-pow(R1,3))/(pow(R2,3)*log(R1)-pow(R1,3)*log(R2));
384-
beta = -3*gammma*(log(R2)-log(R1))/(pow(R1,3)*log(R2)-pow(R2,3)*log(R1)) ;
385-
rho = -(alpha/pow(r,4)*(8*log(r)-6) + 8./3.*beta/r+8*gammma/pow(r,4))*cos(theta) + rho_0;
383+
alpha = -gammma*(std::pow(R2,3)-std::pow(R1,3))/(std::pow(R2,3)*std::log(R1)-std::pow(R1,3)*std::log(R2));
384+
beta = -3*gammma*(std::log(R2)-std::log(R1))/(std::pow(R1,3)*std::log(R2)-std::pow(R2,3)*std::log(R1)) ;
385+
rho = -(alpha/std::pow(r,4)*(8*std::log(r)-6) + 8./3.*beta/r+8*gammma/std::pow(r,4))*std::cos(theta) + rho_0;
386386
}
387387
else
388388
{
389-
alpha=gammma*(mmm+1)*(pow(R1,-3)-pow(R2,-3))/(pow(R1,-mmm-4)-pow(R2,-mmm-4));
390-
beta=-3*gammma*(pow(R1,mmm+1)-pow(R2,mmm+1))/(pow(R1,mmm+4)-pow(R2,mmm+4));
391-
rho= -(2*alpha*pow(r,-4)*(mmm+3)/(mmm+1)*(mmm-1)-2*beta/3*(mmm-1)*(mmm+3)*pow(r,mmm)-mmm*(mmm+5)*2*gammma*pow(r,mmm-3) )*cos(theta) + rho_0;
389+
alpha=gammma*(mmm+1)*(std::pow(R1,-3)-std::pow(R2,-3))/(std::pow(R1,-mmm-4)-std::pow(R2,-mmm-4));
390+
beta=-3*gammma*(std::pow(R1,mmm+1)-std::pow(R2,mmm+1))/(std::pow(R1,mmm+4)-std::pow(R2,mmm+4));
391+
rho= -(2*alpha*std::pow(r,-4)*(mmm+3)/(mmm+1)*(mmm-1)-2*beta/3*(mmm-1)*(mmm+3)*std::pow(r,mmm)-mmm*(mmm+5)*2*gammma*std::pow(r,mmm-3) )*std::cos(theta) + rho_0;
392392
}
393393

394394
out.densities[i] = rho;

benchmarks/king2dcompressible/plugin/code.cc

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -72,13 +72,13 @@ namespace aspect
7272

7373
const double depth = 1.0-position(dim-1);
7474

75-
out.viscosities[i] = ((Di==0.0)?1.0:Di)/Ra*exp(-b*(temperature- this->get_adiabatic_conditions().temperature(position))+c*depth);
75+
out.viscosities[i] = ((Di==0.0)?1.0:Di)/Ra*std::exp(-b*(temperature- this->get_adiabatic_conditions().temperature(position))+c*depth);
7676

7777
out.specific_heat[i] = reference_specific_heat;
7878
out.thermal_conductivities[i] = 1.0;
7979
out.thermal_expansion_coefficients[i] = (Di==0.0)?1.0:Di;
8080

81-
double rho = reference_rho * exp(depth * Di/gamma);
81+
double rho = reference_rho * std::exp(depth * Di/gamma);
8282
rho *= 1.0 - out.thermal_expansion_coefficients[i] * (temperature - this->get_adiabatic_conditions().temperature(position))
8383
+ (tala?0.0:1.0)*Di*gamma
8484
* (pressure - this->get_adiabatic_conditions().pressure(position));

benchmarks/layeredflow/layeredflow.cc

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -50,8 +50,8 @@ namespace aspect
5050
const double yplus=(1.+y0)/beta;
5151
const double yminus=(1.-y0)/beta;
5252
const double C1 = 2.*numbers::PI /
53-
( beta*std::log(beta*beta+pow(1.+y0,2.0))-2.*(1.+y0)*std::atan(yplus)
54-
- beta*std::log(beta*beta+pow(1.-y0,2.0))+2.*(1.-y0)*std::atan(yminus)
53+
( beta*std::log(beta*beta+std::pow(1.+y0,2.0))-2.*(1.+y0)*std::atan(yplus)
54+
- beta*std::log(beta*beta+std::pow(1.-y0,2.0))+2.*(1.-y0)*std::atan(yminus)
5555
+ 2.*numbers::PI*(1.+2.*epsilon) );
5656

5757
const double C2 = ( beta*std::log( beta*beta+Utilities::fixed_power<2>(1+y0) )- 2.*(1.+y0)*std::atan(yplus)

benchmarks/operator_splitting/advection_reaction/advection_reaction.cc

Lines changed: 6 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -126,7 +126,7 @@ namespace aspect
126126
for (unsigned int q=0; q < in.n_evaluation_points(); ++q)
127127
{
128128
// dC/dt = - z * lambda * C
129-
const double decay_constant = half_life > 0.0 ? log(2.0) / half_life : 0.0;
129+
const double decay_constant = half_life > 0.0 ? std::log(2.0) / half_life : 0.0;
130130
const double z = in.position[q](1);
131131
reaction_out->reaction_rates[q][0] = - decay_constant * z / time_scale * in.composition[q][0];
132132
}
@@ -225,7 +225,7 @@ namespace aspect
225225
for (unsigned int q=0; q<heating_model_outputs.heating_source_terms.size(); ++q)
226226
{
227227
// dC/dt = - z * lambda * C
228-
const double decay_constant = half_life > 0.0 ? log(2.0) / half_life : 0.0;
228+
const double decay_constant = half_life > 0.0 ? std::log(2.0) / half_life : 0.0;
229229
const double z = in.position[q](1);
230230
heating_model_outputs.rates_of_temperature_change[q] = - decay_constant * z / time_scale * in.composition[q][0];
231231

@@ -290,8 +290,8 @@ namespace aspect
290290
values[0] = 0.0; // velocity x
291291
values[1] = 0.0; // velocity z
292292
values[2] = 0.0; // pressure
293-
values[3] = sin(2*numbers::PI*(x-t*0.01))*exp(-log(2.0)/10.0*z*t); // temperature
294-
values[4] = sin(2*numbers::PI*(x-t*0.01))*exp(-log(2.0)/10.0*z*t); // composition
293+
values[3] = std::sin(2*numbers::PI*(x-t*0.01))*std::exp(-std::log(2.0)/10.0*z*t); // temperature
294+
values[4] = std::sin(2*numbers::PI*(x-t*0.01))*std::exp(-std::log(2.0)/10.0*z*t); // composition
295295
}
296296
};
297297

@@ -324,7 +324,7 @@ namespace aspect
324324

325325
double t = this->get_time() / time_scale;
326326

327-
return sin(2*numbers::PI*(x-t*0.01))*exp(-log(2.0)/10.0*z*t);
327+
return std::sin(2*numbers::PI*(x-t*0.01))*std::exp(-std::log(2.0)/10.0*z*t);
328328
}
329329

330330

@@ -358,7 +358,7 @@ namespace aspect
358358

359359
double t = this->get_time() / time_scale;
360360

361-
return sin(2*numbers::PI*(x-t*0.01))*exp(-log(2.0)/10.0*z*t);
361+
return std::sin(2*numbers::PI*(x-t*0.01))*std::exp(-std::log(2.0)/10.0*z*t);
362362
}
363363

364364
template <int dim>

benchmarks/operator_splitting/exponential_decay/exponential_decay.cc

Lines changed: 4 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -126,7 +126,7 @@ namespace aspect
126126
for (unsigned int q=0; q < in.n_evaluation_points(); ++q)
127127
{
128128
// dC/dt = - lambda * C
129-
const double decay_constant = half_life > 0.0 ? log(2.0) / half_life : 0.0;
129+
const double decay_constant = half_life > 0.0 ? std::log(2.0) / half_life : 0.0;
130130
reaction_out->reaction_rates[q][0] = - decay_constant / time_scale * in.composition[q][0];
131131
}
132132
}
@@ -225,7 +225,7 @@ namespace aspect
225225
for (unsigned int q=0; q<heating_model_outputs.heating_source_terms.size(); ++q)
226226
{
227227
// dC/dt = - lambda * C
228-
const double decay_constant = half_life > 0.0 ? log(2.0) / half_life : 0.0;
228+
const double decay_constant = half_life > 0.0 ? std::log(2.0) / half_life : 0.0;
229229
heating_model_outputs.rates_of_temperature_change[q] = - decay_constant / time_scale * in.composition[q][0];
230230

231231
heating_model_outputs.heating_source_terms[q] = 0.0;
@@ -285,8 +285,8 @@ namespace aspect
285285
values[0] = 0.0; // velocity x
286286
values[1] = 0.0; // velocity z
287287
values[2] = 0.0; // pressure
288-
values[3] = exp(-log(2.0)/10.0*this->get_time()); // temperature
289-
values[4] = exp(-log(2.0)/10.0*this->get_time()); // composition
288+
values[3] = std::exp(-std::log(2.0)/10.0*this->get_time()); // temperature
289+
values[4] = std::exp(-std::log(2.0)/10.0*this->get_time()); // composition
290290
}
291291
};
292292

benchmarks/shear_bands/shear_bands.cc

Lines changed: 5 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -86,7 +86,7 @@ namespace aspect
8686

8787
double reference_darcy_coefficient () const override
8888
{
89-
return reference_permeability * pow(0.01, permeability_exponent) / eta_f;
89+
return reference_permeability * std::pow(0.01, permeability_exponent) / eta_f;
9090
}
9191

9292

@@ -150,7 +150,7 @@ namespace aspect
150150
{
151151
double porosity = std::max(in.composition[i][porosity_idx],1e-4);
152152

153-
melt_out->compaction_viscosities[i] = xi_0 * pow(porosity/background_porosity,-compaction_viscosity_exponent);
153+
melt_out->compaction_viscosities[i] = xi_0 * std::pow(porosity/background_porosity,-compaction_viscosity_exponent);
154154
melt_out->fluid_viscosities[i]= eta_f;
155155
melt_out->permeabilities[i]= reference_permeability * std::pow(porosity,permeability_exponent);
156156
melt_out->fluid_densities[i]= reference_rho_f;
@@ -534,8 +534,8 @@ namespace aspect
534534
PlaneWaveMeltBandsInitialCondition<dim>::
535535
initial_composition (const Point<dim> &position, const unsigned int /*n_comp*/) const
536536
{
537-
return background_porosity * (1.0 + amplitude * cos(wave_number*position[0]*sin(initial_band_angle)
538-
+ wave_number*position[1]*cos(initial_band_angle)));
537+
return background_porosity * (1.0 + amplitude * std::cos(wave_number*position[0]*std::sin(initial_band_angle)
538+
+ wave_number*position[1]*std::cos(initial_band_angle)));
539539
}
540540

541541

@@ -736,7 +736,7 @@ namespace aspect
736736
const double min_velocity = bm.boundary_velocity(lower_boundary, lower_boundary_point).norm();
737737

738738
const double strain_rate = 0.5 * (max_velocity + min_velocity) / this->get_geometry_model().maximal_depth();
739-
const double theta = std::atan(std::sin(initial_band_angle) / (std::cos(initial_band_angle) - time * strain_rate/sqrt(2.0) * std::sin(initial_band_angle)));
739+
const double theta = std::atan(std::sin(initial_band_angle) / (std::cos(initial_band_angle) - time * strain_rate/numbers::SQRT2 * std::sin(initial_band_angle)));
740740
const double analytical_growth_rate = - eta_0 / (xi_0 + 4.0/3.0 * eta_0) * alpha * (1.0 - background_porosity)
741741
* 2.0 * strain_rate * std::sin(2.0 * theta);
742742

0 commit comments

Comments
 (0)