update 2D implementation

src/convection/incflo_compute_MAC_projected_velocities.cpp

@@ -165,7 +165,7 @@ incflo::compute_MAC_projected_velocities (
         // Predict normal velocity to faces -- note that the {u_mac, v_mac, w_mac}
         //    returned from this call are on face CENTROIDS
         bool allow_inflow_on_outflow = false;
-        /*HydroUtils::ExtrapVelToFaces(*vel[lev], *vel_forces[lev],
+        HydroUtils::ExtrapVelToFaces(*vel[lev], *vel_forces[lev],
                                       AMREX_D_DECL(*u_mac[lev], *v_mac[lev], *w_mac[lev]),
                                       get_velocity_bcrec(), get_velocity_bcrec_device_ptr(),
                                       geom[lev], l_dt,
@@ -175,7 +175,7 @@ incflo::compute_MAC_projected_velocities (
 #endif
                                       m_godunov_ppm, m_godunov_use_forces_in_trans,
                                       l_advection_type, PPM::default_limiter,
-                                      allow_inflow_on_outflow, BC_MF.get());*/
+                                      allow_inflow_on_outflow, BC_MF.get());
 
         //add surface tension
         //if(m_vof_advect_tracer)

src/vof/VolumeOfFluid.cpp

@@ -2566,7 +2566,7 @@ VolumeOfFluid::tracer_vof_init_fraction (int lev, MultiFab& a_tracer)
             Array<Real,AMREX_SPACEDIM> center{AMREX_D_DECL(0.5*(problo[0]+probhi[0]),
                                                            0.5*(problo[1]+probhi[1]),
                                                            0.5*(problo[2]+probhi[2]))};
-            Real radius = .3; //5.0*dx[0];
+            Real radius = .2; //5.0*dx[0];
             bool fluid_is_inside = true;
 
             EB2::SphereIF my_sphere(radius, center, fluid_is_inside);
@@ -2602,7 +2602,7 @@ VolumeOfFluid::tracer_vof_init_fraction (int lev, MultiFab& a_tracer)
     //auto my_box=  EB2::rotate(EB2::BoxIF( low,  high, fluid_is_inside), .3, 1);
     auto my_box1=  EB2::rotate(my_box, .3, 0);
     auto my_box2=  EB2::rotate(my_box1, .2, 2);
-    auto two =EB2::makeIntersection(my_sphere, my_sphere1);
+    auto two =EB2::makeIntersection(my_sphere, my_box);
    //auto two = EB2::makeComplement(EB2::makeUnion(my_cyl_1, my_cyl));
 
     // Generate GeometryShop
@@ -2657,15 +2657,14 @@ VolumeOfFluid::tracer_vof_init_fraction (int lev, MultiFab& a_tracer)
     // Read the file line by line
     std::string line;
     while (std::getline(infile, line)) {
-        std::istringstream iss(line);
-        int i, j, k;
-        Real value;
-        if (!(iss >> i >> j >> k >> value)) {
-            std::cerr << "Error reading line: " << line << std::endl;
-            continue;
-        }
+       std::istringstream iss(line);
+       int i, j, k;
+       Real value;
+       if (!(iss >> i >> j >> k >> value)) {
+          std::cerr << "Error reading line: " << line << std::endl;
+          continue;
+       }
        vout.emplace_back(value,i,j,k);
-
     }
     infile.close();
 #ifdef AMRE_USE_OMP
@@ -2678,12 +2677,21 @@ VolumeOfFluid::tracer_vof_init_fraction (int lev, MultiFab& a_tracer)
             auto const& tracer = a_tracer.array(mfi);
 
             for(int n=0;n<vout.size();++n)
-            /* if(vout[n].i>=vbx.smallEnd()[0]&&vout[n].i<=vbx.bigEnd()[0]&&
-                vout[n].j>=vbx.smallEnd()[1]&&vout[n].j<=vbx.bigEnd()[1]&&
-                vout[n].k>=vbx.smallEnd()[2]&&vout[n].k<=vbx.bigEnd()[2]){*/
-             if(vbx.contains(vout[n].i,vout[n].j,vout[n].k)){
+              if(vout[n].i>=vbx.smallEnd()[0]&&vout[n].i<=vbx.bigEnd()[0]&&
+                 vout[n].j>=vbx.smallEnd()[1]&&vout[n].j<=vbx.bigEnd()[1]
+#if AMREX_SPACEDIM==2
+                 &&vout[n].k==7)
+#else
+                 &&vout[n].k>=vbx.smallEnd()[2]&&vout[n].k<=vbx.bigEnd()[2])
+#endif
+              {
+             //if(vbx.contains(vout[n].i,vout[n].j,vout[n].k)){
+#if AMREX_SPACEDIM==2
+                tracer(vout[n].i,vout[n].j,0,0)=vout[n].vof;
+#else
                 tracer(vout[n].i,vout[n].j,vout[n].k,0)=vout[n].vof;
-             }
+#endif
+              }
 
             /*amrex::ParallelFor(vbx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
             {

test_2d/inputs.droplet

@@ -1,8 +1,8 @@
 #¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
 #            SIMULATION STOP            #
 #.......................................#
-stop_time               =   8.1         # Max (simulated) time to evolve
-max_step                = 20          # Max number of time steps
+stop_time               =   8         # Max (simulated) time to evolve
+max_step                = 10000          # Max number of time steps
 steady_state            =   0           # Steady-state solver? 
 
 #¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
@@ -14,7 +14,7 @@ incflo.cfl              =   0.1        # CFL factor
 #¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
 #            INPUT AND OUTPUT           #
 #.......................................#
-amr.plot_int            =  2           # Steps between plot files
+amr.plot_int            =  5           # Steps between plot files
 amr.check_int           =   1000        # Steps between checkpoint files
 amr.restart             =   ""          # Checkpoint to restart from 
 
@@ -29,12 +29,12 @@ incflo.mu               =   1.e-3         # Dynamic viscosity coefficient
 incflo.vof_advect_tracer=   1
 incflo.mu_s             =   1.0 
 incflo.ro_s             =   1000.0
-incflo.sigma            =   0.
+incflo.sigma            =   10.
 
 #¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
 #        ADAPTIVE MESH REFINEMENT       #
 #.......................................#
-amr.n_cell              =   16 16      # Grid cells at coarsest AMRlevel
+amr.n_cell              =   32 32      # Grid cells at coarsest AMRlevel
 amr.max_level           =   0           # Max AMR level in hierarchy 
 amr.max_grid_size       =   64 64 
 
@@ -43,8 +43,8 @@ amr.max_grid_size       =   64 64
 #.......................................#
 geometry.prob_lo        =   0.  0.     # Lo corner coordinates
 geometry.prob_hi        =   1.  1.    # Hi corner coordinates
-geometry.is_periodic    =   1   1       # Periodicity x y z (0/1)
-
+geometry.is_periodic    =   0 0       # Periodicity x y z (0/1)
+geometry.coord_sys   =  1
 
 #¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
 #           INITIAL CONDITIONS          #
@@ -54,14 +54,14 @@ incflo.ic_u             =   0.
 incflo.ic_v             =   0.
 incflo.ic_w             =   0.
 # Boundary conditions
-#ylo.type                =   "po"
-#ylo.pressure            =   0.  
-#yhi.type                =   "po"
-#yhi.pressure            =   0.  
-#xlo.type                =   "po"
-#xlo.pressure            =   0.
-#xhi.type                =   "po"
-#xhi.pressure            =   0.0
+ylo.type                =   "po"
+ylo.pressure            =   0.  
+yhi.type                =   "po"
+yhi.pressure            =   0.  
+xlo.type                =   "po"
+xlo.pressure            =   0.
+xhi.type                =   "po"
+xhi.pressure            =   0.0
 #zlo.type                =   "po"
 #zlo.pressure            =   0.
 #zhi.type                =   "po"

test_3d/inputs.droplet

@@ -2,7 +2,7 @@
 #            SIMULATION STOP            #
 #.......................................#
 stop_time               =   6         # Max (simulated) time to evolve
-max_step                =  2000          # Max number of time steps
+max_step                =  20000          # Max number of time steps
 steady_state            =   0           # Steady-state solver? 
 
 #¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#