Merge branch 'development' into new_shock_flag_call

Diagnostics/DustCollapse/main.cpp

@@ -102,18 +102,18 @@ int main(int argc, char* argv[])
 #if (AMREX_SPACEDIM == 1)
         int nbins = domain.length(0);
 #elif (AMREX_SPACEDIM == 2)
-        auto x_maxdist = fmax(fabs(problo[0]), fabs(probhi[0]));
-        auto y_maxdist = fmax(fabs(problo[1] - yctr), fabs(probhi[1] - yctr));
+        auto x_maxdist = std::max(std::abs(problo[0]), std::abs(probhi[0]));
+        auto y_maxdist = std::max(std::abs(problo[1] - yctr), std::abs(probhi[1] - yctr));
 
-        auto max_dist = sqrt(x_maxdist*x_maxdist + y_maxdist*y_maxdist);
+        auto max_dist = std::sqrt(x_maxdist*x_maxdist + y_maxdist*y_maxdist);
 
         int nbins = int(max_dist / dx_fine);
 #else
-        auto x_maxdist = fmax(fabs(problo[0] - xctr), fabs(probhi[0] - xctr));
-        auto y_maxdist = fmax(fabs(problo[1] - yctr), fabs(probhi[1] - yctr));
-        auto z_maxdist = fmax(fabs(problo[2] - zctr), fabs(probhi[2] - zctr));
+        auto x_maxdist = std::max(std::abs(problo[0] - xctr), std::abs(probhi[0] - xctr));
+        auto y_maxdist = std::max(std::abs(problo[1] - yctr), std::abs(probhi[1] - yctr));
+        auto z_maxdist = std::max(std::abs(problo[2] - zctr), std::abs(probhi[2] - zctr));
 
-        auto max_dist = sqrt(x_maxdist*x_maxdist + y_maxdist*y_maxdist +
+        auto max_dist = std::sqrt(x_maxdist*x_maxdist + y_maxdist*y_maxdist +
                              z_maxdist*z_maxdist);
 
         int nbins = int(max_dist / dx_fine);
@@ -147,7 +147,7 @@ int main(int argc, char* argv[])
         // over levels, we will compare to the finest level index space to
         // determine if we've already output here
         int mask_size = domain.length().max();
-        Vector<int> imask(pow(mask_size, AMREX_SPACEDIM), 1);
+        Vector<int> imask(std::pow(mask_size, AMREX_SPACEDIM), 1);
 
         // counter
         int cnt = 0;
@@ -210,21 +210,21 @@ int main(int argc, char* argv[])
         // analytic expression for the radius as a function of time t = 0.00
         Real max_dens = 1.e9;
 
-        if (fabs(data.Time()) <= 1.e-8)
+        if (std::abs(data.Time()) <= 1.e-8)
             max_dens = 1.e9;
-        else if (fabs(data.Time() - 0.01) <= 1.e-8)
+        else if (std::abs(data.Time() - 0.01) <= 1.e-8)
             max_dens = 1.043345e9;
-        else if (fabs(data.Time() - 0.02) <= 1.e-8)
+        else if (std::abs(data.Time() - 0.02) <= 1.e-8)
             max_dens = 1.192524e9;
-        else if (fabs(data.Time() - 0.03) <= 1.e-8)
+        else if (std::abs(data.Time() - 0.03) <= 1.e-8)
             max_dens = 1.527201e9;
-        else if (fabs(data.Time() - 0.04) <= 1.e-8)
+        else if (std::abs(data.Time() - 0.04) <= 1.e-8)
             max_dens = 2.312884e9;
-        else if (fabs(data.Time() - 0.05) <= 1.e-8)
+        else if (std::abs(data.Time() - 0.05) <= 1.e-8)
             max_dens = 4.779133e9;
-        else if (fabs(data.Time() - 0.06) <= 1.e-8)
+        else if (std::abs(data.Time() - 0.06) <= 1.e-8)
             max_dens = 24.472425e9;
-        else if (fabs(data.Time() - 0.065) <= 1.e-8)
+        else if (std::abs(data.Time() - 0.065) <= 1.e-8)
             max_dens = 423.447291e9;
         else {
             Print() << "Dont know the maximum density at this time: " << data.Time() <<std::endl;

Diagnostics/Radiation/Radiation_utils.H

@@ -126,7 +126,7 @@ void WriteSlicefile(const int nbins, const Vector<Real> r,
         slicefile << std::setw(w) << r[i];
 
         for (auto it=vars.begin(); it!=vars.end(); ++it) {
-            if(fabs((*it)[i]) < SMALL) (*it)[i] = 0.0;
+            if(std::abs((*it)[i]) < SMALL) (*it)[i] = 0.0;
             slicefile << std::setw(w) << (*it)[i];
         }
 

Diagnostics/Radiation/gaussian_pulse.cpp

@@ -64,9 +64,9 @@ int main(int argc, char* argv[])
 
     // compute the size of the radially-binned array -- we'll do it to
     // the furtherest corner of the domain
-    double x_maxdist = max(fabs(probhi[0] - xctr), fabs(problo[0] - xctr));
-    double y_maxdist = max(fabs(probhi[1] - yctr), fabs(problo[1] - yctr));
-    double maxdist = sqrt(x_maxdist*x_maxdist + y_maxdist*y_maxdist);
+    double x_maxdist = std::max(std::abs(probhi[0] - xctr), std::abs(problo[0] - xctr));
+    double y_maxdist = std::max(std::abs(probhi[1] - yctr), std::abs(problo[1] - yctr));
+    double maxdist = std::sqrt(x_maxdist*x_maxdist + y_maxdist*y_maxdist);
 
     double dx_fine = *(std::min_element(dx.begin(), dx.end()));
 
@@ -102,7 +102,7 @@ int main(int argc, char* argv[])
     // over levels, we will compare to the finest level index space to
     // determine if we've already output here
     int mask_size = domain.length().max();
-    Vector<int> imask(pow(mask_size, AMREX_SPACEDIM), 1);
+    Vector<int> imask(std::pow(mask_size, AMREX_SPACEDIM), 1);
 
     // loop over the data, starting at the finest grid, and if we haven't
     // already stored data in that grid location (according to imask),

Diagnostics/Radiation/lgt_frnt1d.cpp

@@ -53,7 +53,7 @@ int main(int argc, char* argv[])
 
     // compute the size of the radially-binned array -- we'll do it to
     // the furtherest corner of the domain
-    double maxdist = fabs(probhi[0] - problo[0]);
+    double maxdist = std::abs(probhi[0] - problo[0]);
     double dx_fine = *(std::min_element(dx.begin(), dx.end()));
     int nbins = int(maxdist / dx_fine);
 
@@ -93,7 +93,7 @@ int main(int argc, char* argv[])
     // over levels, we will compare to the finest level index space to
     // determine if we've already output here
     int mask_size = domain.length().max();
-    Vector<int> imask(pow(mask_size, AMREX_SPACEDIM), 1);
+    Vector<int> imask(std::pow(mask_size, AMREX_SPACEDIM), 1);
 
     // loop over the data, starting at the finest grid, and if we haven't
     // already stored data in that grid location (according to imask),

Diagnostics/Radiation/rad_shock.cpp

@@ -106,7 +106,7 @@ int main(int argc, char* argv[])
     // over levels, we will compare to the finest level index space to
     // determine if we've already output here
     int mask_size = domain.length().max();
-    Vector<int> imask(pow(mask_size, AMREX_SPACEDIM), 1);
+    Vector<int> imask(std::pow(mask_size, AMREX_SPACEDIM), 1);
 
     // counter
     auto cnt = 0;
@@ -214,7 +214,7 @@ int main(int argc, char* argv[])
             slicefile << std::setw(w) << vars_bin[isv[i]+((*it)+1)*nbins];
 
         auto it = comps.end()-1;
-        slicefile << std::setw(w) << pow(vars_bin[isv[i]+((*it)+1)*nbins] / arad, 0.25);
+        slicefile << std::setw(w) << std::pow(vars_bin[isv[i]+((*it)+1)*nbins] / arad, 0.25);
 
         slicefile << std::endl;
     }

Diagnostics/Radiation/rad_sphere.cpp

@@ -142,7 +142,7 @@ int main(int argc, char* argv[])
     // over levels, we will compare to the finest level index space to
     // determine if we've already output here
     int mask_size = domain.length()[0];
-    Vector<int> imask(pow(mask_size, AMREX_SPACEDIM), 1);
+    Vector<int> imask(std::pow(mask_size, AMREX_SPACEDIM), 1);
 
     // counter
     int cnt = 0;

Diagnostics/Sedov/main.cpp

@@ -148,7 +148,7 @@ int main(int argc, char* argv[])
 
 #else
     double x_maxdist = amrex::max(std::abs(probhi[0] - xctr),
-                                  std::fabs(problo[0] - xctr));
+                                  std::abs(problo[0] - xctr));
     double y_maxdist = amrex::max(std::abs(probhi[1] - yctr),
                                   std::abs(problo[1] - yctr));
     double z_maxdist = amrex::max(std::abs(probhi[2] - zctr),
@@ -342,10 +342,10 @@ int main(int argc, char* argv[])
     // write the data in columns
     const auto SMALL = 1.e-20;
     for (auto i = 0; i < nbins; i++) {
-        if (fabs(dens_bin[i]) < SMALL) dens_bin[i] = 0.0;
-        if (fabs( vel_bin[i]) < SMALL) vel_bin[i] = 0.0;
-        if (fabs(pres_bin[i]) < SMALL) pres_bin[i] = 0.0;
-        if (fabs(   e_bin[i]) < SMALL) e_bin[i] = 0.0;
+        if (std::abs(dens_bin[i]) < SMALL) dens_bin[i] = 0.0;
+        if (std::abs( vel_bin[i]) < SMALL) vel_bin[i] = 0.0;
+        if (std::abs(pres_bin[i]) < SMALL) pres_bin[i] = 0.0;
+        if (std::abs(   e_bin[i]) < SMALL) e_bin[i] = 0.0;
 
         slicefile << std::setw(w) << r[i] << std::setw(w) << dens_bin[i] << std::setw(w) << vel_bin[i] << std::setw(w) << pres_bin[i] << std::setw(w) << e_bin[i] << std::endl;
     }

Exec/gravity_tests/DustCollapse/problem_bc_fill.H

@@ -5,8 +5,8 @@ AMREX_GPU_HOST_DEVICE AMREX_INLINE
 void problem_bc_fill(int i, int j, int k,
                      Array4<Real> const& state,
                      Real time,
-                     Array1D<BCRec, 0, NUM_STATE-1> bcs,
-                     GeometryData geomdata)
+                     const Array1D<BCRec, 0, NUM_STATE-1>& bcs,
+                     const GeometryData& geomdata)
 {
 #if AMREX_SPACEDIM == 3
     const Real* dx = geomdata.CellSize();

Exec/gravity_tests/hydrostatic_adjust/problem_bc_fill.H

@@ -5,8 +5,8 @@ AMREX_GPU_HOST_DEVICE AMREX_INLINE
 void problem_bc_fill(int i, int j, int k,
                      Array4<Real> const& state,
                      Real time,
-                     Array1D<BCRec, 0, NUM_STATE-1> bcs,
-                     GeometryData geomdata)
+                     const Array1D<BCRec, 0, NUM_STATE-1>& bcs,
+                     const GeometryData& geomdata)
 {
     // Set the velocity to zero at the top physical boundary.
 

Exec/hydro_tests/Noh/problem_bc_fill.H

@@ -8,8 +8,8 @@ AMREX_GPU_HOST_DEVICE AMREX_INLINE
 void problem_bc_fill(int i, int j, int k,
                      Array4<Real> const& state,
                      Real time,
-                     Array1D<BCRec, 0, NUM_STATE-1> bcs,
-                     GeometryData geomdata)
+                     const Array1D<BCRec, 0, NUM_STATE-1>& bcs,
+                     const GeometryData& geomdata)
 {
     const Real pres_init = 1.0e-6_rt;
     const Real rho_init = 1.0e0_rt;

Exec/hydro_tests/double_mach_reflection/problem_bc_fill.H

@@ -5,8 +5,8 @@ AMREX_GPU_HOST_DEVICE AMREX_INLINE
 void problem_bc_fill(int i, int j, int k,
                      Array4<Real> const& state,
                      Real time,
-                     Array1D<BCRec, 0, NUM_STATE-1> bcs,
-                     GeometryData geomdata)
+                     const Array1D<BCRec, 0, NUM_STATE-1>& bcs,
+                     const GeometryData& geomdata)
 {
     const int* domlo = geomdata.Domain().loVect();
     const int* domhi = geomdata.Domain().hiVect();

Exec/radiation_tests/Rad2Tshock/problem_bc_fill.H

@@ -5,8 +5,8 @@ AMREX_GPU_HOST_DEVICE AMREX_INLINE
 void problem_bc_fill(int i, int j, int k,
                      Array4<Real> const& state,
                      Real time,
-                     Array1D<BCRec, 0, NUM_STATE-1> bcs,
-                     GeometryData geomdata)
+                     const Array1D<BCRec, 0, NUM_STATE-1>& bcs,
+                     const GeometryData& geomdata)
 {
     // The strategy here is to set Dirichlet condition for inflow and
     // outflow boundaries, and let the Riemann solver sort out the

Exec/science/Detonation/problem_bc_fill.H

@@ -45,8 +45,8 @@ AMREX_GPU_HOST_DEVICE AMREX_INLINE
 void problem_bc_fill(int i, int j, int k,
                      Array4<Real> const& state,
                      Real time,
-                     Array1D<BCRec, 0, NUM_STATE-1> bcs,
-                     GeometryData geomdata)
+                     const Array1D<BCRec, 0, NUM_STATE-1>& bcs,
+                     const GeometryData& geomdata)
 {
     // Override the generic routine at the physical boundaries by
     // setting the material to the ambient state. Note that we

Exec/science/nova/problem_initialize_state_data.H

@@ -121,10 +121,10 @@ void problem_initialize_state_data (int i, int j, int k,
 
         if (problem::num_vortices % 2 == 0){
             temp_pert = delta * (1.0_rt + std::sin(static_cast<Real>(problem::num_vortices)*M_PI*x/problem::L) )*
-                                    delta*exp(-pow(ydist / problem::width, 2));
+                                    delta*std::exp(-pow(ydist / problem::width, 2));
         } else {
             temp_pert = delta * (1.0_rt + std::cos(static_cast<Real>(problem::num_vortices)*M_PI*x/problem::L) )*
-                                     delta*exp(-pow(ydist / problem::width, 2));
+                                     delta*std::exp(-pow(ydist / problem::width, 2));
         }
 
         state(i,j,k,UTEMP)  = state(i,j,k, UTEMP)*(1.0_rt + temp_pert);

Exec/science/planet/HotJupiter.cpp

@@ -19,10 +19,10 @@ long double Temp(long double P1,long double Pc,long double P_char,long double T_
   long double Temp;
 
   if(Pc>P1){
-    Temp = Td * pow(1.0e0 + (Lad / (Lin - Lad))*pow((P1/Pc),exp1),exp2);
+    Temp = Td * std::pow(1.0e0 + (Lad / (Lin - Lad)) * std::pow((P1/Pc),exp1),exp2);
   }
   else{
-    Temp = T_char * pow(P1/P_char,Lad);
+    Temp = T_char * std::pow(P1/P_char,Lad);
   }
   return Temp;
 }
@@ -95,9 +95,9 @@ int main(){
   //Defining the properties of the planet atmosphere
                             // the density limit at the top of the atmosphere
   P_top=den_top*R*Td ;                       // the pressure limit at the density limit
-  Tc = Td * pow(Lin/(Lin-Lad),exp2);                 //the temperature at the radiative-convective boundary (RCB)
-  Pc = P_char * pow(Tc/T_char,1.0/Lad);              //the pressure at the RCB
-  Pd = pow((Lin-Lad)/(2.0*Lin-Lad),1.0/exp1)*Pc;     //the characteristic pressure defining the radiative zone
+  Tc = Td * std::pow(Lin/(Lin-Lad),exp2);                 //the temperature at the radiative-convective boundary (RCB)
+  Pc = P_char * std::pow(Tc/T_char,1.0/Lad);              //the pressure at the RCB
+  Pd = std::pow((Lin-Lad)/(2.0*Lin-Lad),1.0/exp1)*Pc;     //the characteristic pressure defining the radiative zone
 
 
 
@@ -112,7 +112,7 @@ int main(){
   }
 
   buffer_height=z_top*2.25/4.0;               // the height at which the buffer and the top of the atmosphere meet (not continuous at RCB).-> needs sponge
-  optical_depth_buffer=6.35e-3*T_buffer*pow(den_buffer,2.0)*(z_top-buffer_height);// the optical depth in the buffer region [dimensionless]
+  optical_depth_buffer=6.35e-3*T_buffer*std::pow(den_buffer,2.0)*(z_top-buffer_height);// the optical depth in the buffer region [dimensionless]
 
 
 
@@ -157,28 +157,28 @@ int main(){
       den2 = density(P2, T2, R);
       T3 = T1;
       den3 = den1;
-      dz2 = abs ((-P1+P2)* (1.0/grav)/((1.0/2.0)*(den1+den2)));
-      dz3 = abs ((-P1+P3)* (1.0/grav)/((1.0/2.0)*(den1+den3)));
+      dz2 = std::abs ((-P1+P2)* (1.0/grav)/((1.0/2.0)*(den1+den2)));
+      dz3 = std::abs ((-P1+P3)* (1.0/grav)/((1.0/2.0)*(den1+den3)));
 
     a:  if(difference(dz2, delta_z)*difference(dz3, delta_z)<0.0){
         P_sol = (P2+P3)/2.0;
         T_sol = Temp(P_sol, Pc, P_char, T_char, Td, Lin, Lad, exp1, exp2);
         den_sol = density(P_sol,T_sol, R);
-        dz_sol = abs ((-P1+P_sol)* (1.0/grav)/((1.0/2.0)*(den1+den_sol)));
+        dz_sol = std::abs ((-P1+P_sol)* (1.0/grav)/((1.0/2.0)*(den1+den_sol)));
         count2++;
         }
         else{
           if(difference(dz2, delta_z)>0.0){
           P3 = P3 - P3/dP_factor/1.0;
           T3 = Temp(P3, Pc, P_char, T_char, Td, Lin, Lad, exp1, exp2);
           den3 = density(P3, T3, R);
-          dz3 = abs ((-P1+P3)* (1.0/grav)/((1.0/2.0)*(den1+den3)));
+          dz3 = std::abs ((-P1+P3)* (1.0/grav)/((1.0/2.0)*(den1+den3)));
           }
         else if(difference(dz2, delta_z)<0.0){
           P2 = P2 + P2/dP_factor/1.0;
           T2 = Temp(P2, Pc, P_char, T_char, Td, Lin, Lad, exp1, exp2);
           den2 = density(P2,T2,R);
-          dz2 = abs ((-P1+P2)* (1.0/grav)/((1.0/2.0)*(den1+den2)));
+          dz2 = std::abs ((-P1+P2)* (1.0/grav)/((1.0/2.0)*(den1+den2)));
           }
           goto a;
         }
@@ -192,11 +192,11 @@ int main(){
       if(count2>100000){
         goto b;
       }
-      if(abs(difference(dz_sol, delta_z)) >= error_expected){
+      if(std::abs(difference(dz_sol, delta_z)) >= error_expected){
         goto a;
       }
 
-    b:  cout<<count_cell<<' ' <<setprecision(5) << "Final value, "<<' '<<"height[km]="<< ' '<< z/1.e5<<' ' <<"pressure[bar]=" <<' ' << P_sol/1.e6  <<' '<< "error="<<' ' << abs(difference(dz_sol, delta_z))<<' ' << abs ((-P1+P_sol)/dz_sol+ grav*((1.0/2.0)*(den1+den_sol))) << '\n';
+    b:  cout<<count_cell<<' ' <<setprecision(5) << "Final value, "<<' '<<"height[km]="<< ' '<< z/1.e5<<' ' <<"pressure[bar]=" <<' ' << P_sol/1.e6  <<' '<< "error="<<' ' << std::abs(difference(dz_sol, delta_z))<<' ' << std::abs ((-P1+P_sol)/dz_sol+ grav*((1.0/2.0)*(den1+den_sol))) << '\n';
       T[count] = T_sol;
       P[count] = P_sol;
       den[count] = den_sol;
@@ -207,9 +207,9 @@ int main(){
       P1 = P_sol;
       count_cell++;
       if(P_sol < Pc){
-      optical_depth_radiative = optical_depth_radiative+kappa0*pow(den_sol,2.0)*T_sol*delta_z;
+      optical_depth_radiative = optical_depth_radiative+kappa0*std::pow(den_sol,2.0)*T_sol*delta_z;
       }
-      optical_depth_atm = optical_depth_atm+kappa0*pow(den_sol,2.0)*T_sol*delta_z;
+      optical_depth_atm = optical_depth_atm+kappa0*std::pow(den_sol,2.0)*T_sol*delta_z;
       if(optical_depth_atm+optical_depth_buffer>=1.0 && P_at_1opticaldepth==0.0 ){
         P_at_1opticaldepth=P_sol;
         T_at_1opticaldepth=T_sol;

Exec/science/wdmerger/Prob.cpp

@@ -252,7 +252,7 @@ Castro::wd_update (Real time, Real dt)
     bool local_flag = true;
 
     for (int i = 0; i <= 6; ++i) {
-        Castro::volInBoundary(time, vol_P[i], vol_S[i], pow(10.0,i), local_flag);
+        Castro::volInBoundary(time, vol_P[i], vol_S[i], std::pow(10.0,i), local_flag);
     }
 
     // Do all of the reductions.
@@ -355,12 +355,12 @@ Castro::wd_update (Real time, Real dt)
 
     if (mass_P > 0.0 && vol_P[2] > 0.0) {
         rho_avg_P = mass_P / vol_P[2];
-        t_ff_P = sqrt(3.0 * M_PI / (32.0 * C::Gconst * rho_avg_P));
+        t_ff_P = std::sqrt(3.0 * M_PI / (32.0 * C::Gconst * rho_avg_P));
     }
 
     if (mass_S > 0.0 && vol_S[2] > 0.0) {
         rho_avg_S = mass_S / vol_S[2];
-        t_ff_S = sqrt(3.0 * M_PI / (32.0 * C::Gconst * rho_avg_S));
+        t_ff_S = std::sqrt(3.0 * M_PI / (32.0 * C::Gconst * rho_avg_S));
     }
 
     // Compute updated roche Radii

Exec/unit_tests/particles_test/Prob.cpp

@@ -100,7 +100,7 @@ void Castro::problem_post_simulation(Vector<std::unique_ptr<AmrLevel> >& amr_lev
             auto deltax = it->m_pos[0] - p.m_pos[0];
             auto deltay = it->m_pos[1] - p.m_pos[1];
 
-            auto delta = sqrt(deltax*deltax + deltay*deltay);
+            auto delta = std::sqrt(deltax*deltax + deltay*deltay);
 
             // quick nan check here
             if (delta == delta)

Source/driver/Castro.H

@@ -1487,8 +1487,10 @@ protected:
 /// sdc
 ///
     int sdc_iteration;
-    int current_sdc_node;
 
+#ifdef TRUE_SDC
+    int current_sdc_node;
+#endif
 
 
 /* problem-specific includes */