Skip to content

Commit

Permalink
only update corners after doing all boundaries and change logic to be…
Browse files Browse the repository at this point in the history
… if not periodic
  • Loading branch information
hklion committed Aug 30, 2024
1 parent 7bf2aac commit 2d2e8a4
Show file tree
Hide file tree
Showing 7 changed files with 119 additions and 73 deletions.
10 changes: 5 additions & 5 deletions Source/BoundaryConditions/BoundaryConditions_netcdf.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -79,7 +79,6 @@ REMORA::fill_from_bdyfiles (MultiFab& mf_to_fill, const MultiFab& mf_mask, const
const bool apply_south = domain_bcs_type[bccomp+icomp].lo(1) == REMORABCType::clamped;
const bool apply_north = domain_bcs_type[bccomp+icomp].hi(1) == REMORABCType::clamped;


#ifdef AMREX_USE_OMP
#pragma omp parallel if (Gpu::notInLaunchRegion())
#endif
Expand Down Expand Up @@ -139,25 +138,26 @@ REMORA::fill_from_bdyfiles (MultiFab& mf_to_fill, const MultiFab& mf_mask, const
});
}

if (!xlo_ylo.isEmpty() && apply_east && apply_south) {
// If we've applied boundary conditions to either side, update the corner
if (!xlo_ylo.isEmpty() && (apply_east || apply_south)) {
ParallelFor(xlo_ylo, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k,icomp+icomp_to_fill) = 0.5 * (dest_arr(i,dom_lo.y+mf_index_type[1],k,icomp+icomp_to_fill) + dest_arr(dom_lo.x+mf_index_type[0],j,k,icomp+icomp_to_fill));
});
}
if (!xlo_yhi.isEmpty() && apply_east && apply_north) {
if (!xlo_yhi.isEmpty() && (apply_east || apply_north)) {
ParallelFor(xlo_yhi, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k,icomp+icomp_to_fill) = 0.5 * (dest_arr(i,dom_hi.y-mf_index_type[1],k,icomp+icomp_to_fill) + dest_arr(dom_lo.x+mf_index_type[0],j,k,icomp+icomp_to_fill));
});
}
if (!xhi_ylo.isEmpty() && apply_west && apply_south) {
if (!xhi_ylo.isEmpty() && (apply_west || apply_south)) {
ParallelFor(xhi_ylo, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k,icomp+icomp_to_fill) = 0.5 * (dest_arr(i,dom_lo.y+mf_index_type[1],k,icomp+icomp_to_fill) + dest_arr(dom_hi.x-mf_index_type[0],j,k,icomp+icomp_to_fill));
});
}
if (!xhi_yhi.isEmpty() && apply_west && apply_north) {
if (!xhi_yhi.isEmpty() && (apply_west || apply_north)) {
ParallelFor(xhi_yhi, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k,icomp+icomp_to_fill) = 0.5 * (dest_arr(i,dom_hi.y-mf_index_type[1],k,icomp+icomp_to_fill) + dest_arr(dom_hi.x-mf_index_type[0],j,k,icomp+icomp_to_fill));
Expand Down
46 changes: 30 additions & 16 deletions Source/BoundaryConditions/BoundaryConditions_xvel.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -171,22 +171,36 @@ void REMORAPhysBCFunct::impose_xvel_bcs (const Array4<Real>& dest_arr, const Box
Box xlo_yhi = xlo & yhi;
Box xhi_ylo = xhi & ylo;
Box xhi_yhi = xhi & yhi;
ParallelFor(xlo_ylo, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k) = 0.5 * (dest_arr(i,dom_lo.y,k) + dest_arr(dom_lo.x+1,j,k));
});
ParallelFor(xlo_yhi, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k) = 0.5 * (dest_arr(i,dom_hi.y,k) + dest_arr(dom_lo.x+1,j,k));
});
ParallelFor(xhi_ylo, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k) = 0.5 * (dest_arr(i,dom_lo.y,k) + dest_arr(dom_hi.x,j,k));
});
ParallelFor(xhi_yhi, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k) = 0.5 * (dest_arr(i,dom_hi.y,k) + dest_arr(dom_hi.x,j,k));
});

const bool clamp_east = m_domain_bcs_type[bccomp].lo(0) == REMORABCType::clamped;
const bool clamp_west = m_domain_bcs_type[bccomp].hi(0) == REMORABCType::clamped;
const bool clamp_south = m_domain_bcs_type[bccomp].lo(1) == REMORABCType::clamped;
const bool clamp_north = m_domain_bcs_type[bccomp].hi(1) == REMORABCType::clamped;

if (!clamp_east && !clamp_south) {
ParallelFor(xlo_ylo, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k) = 0.5 * (dest_arr(i,dom_lo.y,k) + dest_arr(dom_lo.x+1,j,k));
});
}
if (!clamp_east && !clamp_north) {
ParallelFor(xlo_yhi, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k) = 0.5 * (dest_arr(i,dom_hi.y,k) + dest_arr(dom_lo.x+1,j,k));
});
}
if (!clamp_west && !clamp_south) {
ParallelFor(xhi_ylo, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k) = 0.5 * (dest_arr(i,dom_lo.y,k) + dest_arr(dom_hi.x,j,k));
});
}
if (!clamp_west && !clamp_north) {
ParallelFor(xhi_yhi, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k) = 0.5 * (dest_arr(i,dom_hi.y,k) + dest_arr(dom_hi.x,j,k));
});
}
}


Expand Down
45 changes: 29 additions & 16 deletions Source/BoundaryConditions/BoundaryConditions_yvel.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -173,22 +173,35 @@ void REMORAPhysBCFunct::impose_yvel_bcs (const Array4<Real>& dest_arr, const Box
Box xlo_yhi = xlo & yhi;
Box xhi_ylo = xhi & ylo;
Box xhi_yhi = xhi & yhi;
ParallelFor(xlo_ylo, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k) = 0.5 * (dest_arr(i,dom_lo.y+1,k) + dest_arr(dom_lo.x,j,k));
});
ParallelFor(xlo_yhi, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k) = 0.5 * (dest_arr(i,dom_hi.y,k) + dest_arr(dom_lo.x,j,k));
});
ParallelFor(xhi_ylo, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k) = 0.5 * (dest_arr(i,dom_lo.y+1,k) + dest_arr(dom_hi.x,j,k));
});
ParallelFor(xhi_yhi, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k) = 0.5 * (dest_arr(i,dom_hi.y,k) + dest_arr(dom_hi.x,j,k));
});
const bool clamp_east = m_domain_bcs_type[bccomp].lo(0) == REMORABCType::clamped;
const bool clamp_west = m_domain_bcs_type[bccomp].hi(0) == REMORABCType::clamped;
const bool clamp_south = m_domain_bcs_type[bccomp].lo(1) == REMORABCType::clamped;
const bool clamp_north = m_domain_bcs_type[bccomp].hi(1) == REMORABCType::clamped;

if (!clamp_east && !clamp_south) {
ParallelFor(xlo_ylo, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k) = 0.5 * (dest_arr(i,dom_lo.y+1,k) + dest_arr(dom_lo.x,j,k));
});
}
if (!clamp_east && !clamp_north) {
ParallelFor(xlo_yhi, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k) = 0.5 * (dest_arr(i,dom_hi.y,k) + dest_arr(dom_lo.x,j,k));
});
}
if (!clamp_west && !clamp_south) {
ParallelFor(xhi_ylo, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k) = 0.5 * (dest_arr(i,dom_lo.y+1,k) + dest_arr(dom_hi.x,j,k));
});
}
if (!clamp_west && !clamp_north) {
ParallelFor(xhi_yhi, [=] AMREX_GPU_DEVICE (int i, int j, int k)
{
dest_arr(i,j,k) = 0.5 * (dest_arr(i,dom_hi.y,k) + dest_arr(dom_hi.x,j,k));
});
}
}

Gpu::streamSynchronize();
Expand Down
32 changes: 20 additions & 12 deletions Source/TimeIntegration/REMORA_gls.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -34,6 +34,10 @@ REMORA::gls_prestep (int lev, MultiFab* mf_gls, MultiFab* mf_tke,
const auto dlo = amrex::lbound(domain);
const auto dhi = amrex::ubound(domain);

GeometryData const& geomdata = geom[0].data();
bool is_periodic_in_x = geomdata.isPeriodic(0);
bool is_periodic_in_y = geomdata.isPeriodic(1);

int ncomp = 1;
Vector<BCRec> bcrs_x(ncomp);
Vector<BCRec> bcrs_y(ncomp);
Expand Down Expand Up @@ -77,11 +81,11 @@ REMORA::gls_prestep (int lev, MultiFab* mf_gls, MultiFab* mf_tke,

// Adjust boundaries
// TODO: Make sure indices match with what ROMS does
if (i == dlo.x-1 && (bcr_x.lo(0) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
if (i == dlo.x-1 && !is_periodic_in_x) {
grad_im1 = tke(i,j,k,nstp) - tke(i-1,j,k,nstp);
gradL_im1 = gls(i,j,k,nstp) - gls(i-1,j,k,nstp);
}
else if (i == dhi.x+1 && (bcr_x.hi(0) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
else if (i == dhi.x+1 && !is_periodic_in_x) {
grad_ip1 = tke(i,j,k,nstp) - tke(i-1,j,k,nstp);
gradL_ip1 = gls(i,j,k,nstp) - gls(i-1,j,k,nstp);
}
Expand All @@ -104,11 +108,11 @@ REMORA::gls_prestep (int lev, MultiFab* mf_gls, MultiFab* mf_tke,

// Adjust boundaries
// TODO: Make sure indices match with what ROMS does
if (j == dlo.y-1 && (bcr_y.lo(1) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
if (j == dlo.y-1 && !is_periodic_in_y) {
grad_jm1 = tke(i,j,k,nstp) - tke(i,j-1,k,nstp);
gradL_jm1 = gls(i,j,k,nstp) - gls(i,j-1,k,nstp);
}
else if (j == dhi.y+1 && (bcr_y.hi(1) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
else if (j == dhi.y+1 && !is_periodic_in_y) {
grad_jp1 = tke(i,j,k,nstp) - tke(i,j-1,k,nstp);
gradL_jp1 = gls(i,j,k,nstp) - gls(i,j-1,k,nstp);
}
Expand Down Expand Up @@ -371,6 +375,14 @@ REMORA::gls_corrector (int lev, MultiFab* mf_gls, MultiFab* mf_tke,
MultiFab mf_BCK(ba,dm,1,IntVect(NGROW,NGROW,0));
MultiFab mf_BCP(ba,dm,1,IntVect(NGROW,NGROW,0));

const Box& domain = geom[0].Domain();
const auto dlo = amrex::lbound(domain);
const auto dhi = amrex::ubound(domain);

GeometryData const& geomdata = geom[0].data();
bool is_periodic_in_x = geomdata.isPeriodic(0);
bool is_periodic_in_y = geomdata.isPeriodic(1);


for ( MFIter mfi(*mf_gls, TilingIfNotGPU()); mfi.isValid(); ++mfi )
{
Expand All @@ -391,10 +403,6 @@ REMORA::gls_corrector (int lev, MultiFab* mf_gls, MultiFab* mf_tke,
auto CF = mf_CF.array(mfi);
auto shear2_cached = mf_shear2_cached.array(mfi);

const Box& domain = geom[0].Domain();
const auto dlo = amrex::lbound(domain);
const auto dhi = amrex::ubound(domain);

ParallelFor(gbx1D, [=] AMREX_GPU_DEVICE (int i, int j, int )
{
CF(i,j,0) = 0.0_rt;
Expand Down Expand Up @@ -532,12 +540,12 @@ REMORA::gls_corrector (int lev, MultiFab* mf_gls, MultiFab* mf_tke,
{
Real gradK, gradK_ip1, gradP, gradP_ip1;

if (i == dlo.x-1 && (bcr_x.lo(0) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
if (i == dlo.x-1 && !is_periodic_in_x) {
gradK_ip1 = tke(i+1,j,k,2)-tke(i ,j,k,2);
gradK = gradK_ip1;
gradP_ip1 = gls(i+1,j,k,2)-gls(i ,j,k,2);
gradP = gradP_ip1;
} else if (i == dhi.x+1 && (bcr_x.hi(0) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
} else if (i == dhi.x+1 && !is_periodic_in_x) {
gradK = tke(i ,j,k,2)-tke(i-1,j,k,2);
gradK_ip1 = gradK;
gradP = gls(i ,j,k,2)-gls(i-1,j,k,2);
Expand Down Expand Up @@ -570,11 +578,11 @@ REMORA::gls_corrector (int lev, MultiFab* mf_gls, MultiFab* mf_tke,
Real gradP = (gls(i,j ,k,2)-gls(i,j-1,k,2)) * mskv(i,j ,0);
Real gradP_jp1 = (gls(i,j+1,k,2)-gls(i,j ,k,2)) * mskv(i,j+1,0);

if (j == dlo.y-1 && (bcr_y.lo(1) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
if (j == dlo.y-1 && !is_periodic_in_y) {
gradK = gradK_jp1;
gradP = gradP_jp1;
}
else if (j == dhi.y+1 && (bcr_y.hi(1) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
else if (j == dhi.y+1 && !is_periodic_in_y) {
gradK_jp1 = gradK;
gradP_jp1 = gradP;
}
Expand Down
20 changes: 12 additions & 8 deletions Source/TimeIntegration/REMORA_rhs_uv_2d.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -30,6 +30,10 @@ REMORA::rhs_uv_2d (const Box& xbx, const Box& ybx,
const auto dlo = amrex::lbound(domain);
const auto dhi = amrex::ubound(domain);

GeometryData const& geomdata = geom[0].data();
bool is_periodic_in_x = geomdata.isPeriodic(0);
bool is_periodic_in_y = geomdata.isPeriodic(1);

int ncomp = 1;
Vector<BCRec> bcrs_x(ncomp);
Vector<BCRec> bcrs_y(ncomp);
Expand Down Expand Up @@ -82,11 +86,11 @@ REMORA::rhs_uv_2d (const Box& xbx, const Box& ybx,
Real Huxx_i = DUon(i-1,j,0)-2.0_rt*DUon(i ,j,0)+DUon(i+1,j,0);
Real Huxx_ip1 = DUon(i ,j,0)-2.0_rt*DUon(i+1,j,0)+DUon(i+2,j,0);

if (i == dlo.x && (bcr_x.lo(0) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
if (i == dlo.x && !is_periodic_in_x) {
uxx_i = uxx_ip1;
Huxx_i = Huxx_ip1;
}
else if (i == dhi.x && (bcr_x.hi(0) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
else if (i == dhi.x && !is_periodic_in_x) {
uxx_ip1 = uxx_i;
Huxx_ip1 = Huxx_i;
}
Expand All @@ -110,9 +114,9 @@ REMORA::rhs_uv_2d (const Box& xbx, const Box& ybx,
Real Hvxx_i = DVom(i-1,j,0)-2.0_rt*DVom(i ,j,0)+DVom(i+1,j,0);
Real Hvxx_im1 = DVom(i-2,j,0)-2.0_rt*DVom(i-1,j,0)+DVom(i ,j,0);

if (j == dlo.y and (bcr_y.lo(1) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
if (j == dlo.y and !is_periodic_in_y) {
uee_jm1 = uee_j;
} else if (j == dhi.y+1 and (bcr_y.hi(1) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
} else if (j == dhi.y+1 and !is_periodic_in_y) {
uee_j = uee_jm1;
}

Expand Down Expand Up @@ -174,9 +178,9 @@ REMORA::rhs_uv_2d (const Box& xbx, const Box& ybx,
Real Huee_j = DUon(i,j-1,0)-2.0_rt*DUon(i,j ,0)+DUon(i,j+1,0);
Real Huee_jm1 = DUon(i,j-2,0)-2.0_rt*DUon(i,j-1,0)+DUon(i,j ,0);

if (i == dlo.x and (bcr_x.lo(0) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
if (i == dlo.x and !is_periodic_in_x) {
vxx_im1 = vxx_i;
} else if (i == dhi.x + 1 and (bcr_x.hi(0) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
} else if (i == dhi.x + 1 and !is_periodic_in_x) {
vxx_i = vxx_im1;
}

Expand All @@ -195,11 +199,11 @@ REMORA::rhs_uv_2d (const Box& xbx, const Box& ybx,
Real Hvee_j = DVom(i,j-1,0)-2.0_rt*DVom(i,j ,0)+DVom(i,j+1,0);
Real Hvee_jp1 = DVom(i,j ,0)-2.0_rt*DVom(i,j+1,0)+DVom(i,j+2,0);

if (j == dlo.y and (bcr_y.lo(1) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
if (j == dlo.y and !is_periodic_in_y) {
vee_j = vee_jp1;
Hvee_j = Hvee_jp1;
}
else if (j == dhi.y and (bcr_y.hi(1) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
else if (j == dhi.y and !is_periodic_in_y) {
vee_jp1 = vee_j;
Hvee_jp1 = Hvee_j;
}
Expand Down
20 changes: 12 additions & 8 deletions Source/TimeIntegration/REMORA_rhs_uv_3d.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -51,6 +51,10 @@ REMORA::rhs_uv_3d (const Box& xbx, const Box& ybx,
const auto dlo = amrex::lbound(domain);
const auto dhi = amrex::ubound(domain);

GeometryData const& geomdata = geom[0].data();
bool is_periodic_in_x = geomdata.isPeriodic(0);
bool is_periodic_in_y = geomdata.isPeriodic(1);

int ncomp = 1;
Vector<BCRec> bcrs_x(ncomp);
Vector<BCRec> bcrs_y(ncomp);
Expand Down Expand Up @@ -108,11 +112,11 @@ REMORA::rhs_uv_3d (const Box& xbx, const Box& ybx,
Real Huxx_i = Huon(i-1,j,k)-2.0_rt*Huon(i ,j,k)+Huon(i+1,j,k);
Real Huxx_ip1 = Huon(i ,j,k)-2.0_rt*Huon(i+1,j,k)+Huon(i+2,j,k);

if (i == dlo.x && (bcr_x.lo(0) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
if (i == dlo.x && !is_periodic_in_x) {
uxx_i = uxx_ip1;
Huxx_i = Huxx_ip1;
}
else if (i == dhi.x && (bcr_x.hi(0) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
else if (i == dhi.x && !is_periodic_in_x) {
uxx_ip1 = uxx_i;
Huxx_ip1 = Huxx_i;
}
Expand All @@ -135,9 +139,9 @@ REMORA::rhs_uv_3d (const Box& xbx, const Box& ybx,
Real uee_jm1 = uold(i,j-2,k,nrhs) - 2.0_rt*uold(i,j-1,k,nrhs) + uold(i ,j,k,nrhs);
Real uee_j = uold(i,j-1,k,nrhs) - 2.0_rt*uold(i,j ,k,nrhs) + uold(i,j+1,k,nrhs);

if (j == dlo.y and (bcr_y.lo(1) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
if (j == dlo.y and !is_periodic_in_y) {
uee_jm1 = uee_j;
} else if (j == dhi.y+1 and (bcr_y.hi(1) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
} else if (j == dhi.y+1 and !is_periodic_in_y) {
uee_j = uee_jm1;
}

Expand Down Expand Up @@ -263,9 +267,9 @@ REMORA::rhs_uv_3d (const Box& xbx, const Box& ybx,
Real vxx_im1 = vold(i-2,j,k,nrhs)-2.0_rt*vold(i-1,j,k,nrhs)+vold(i ,j,k,nrhs);
Real vxx_i = vold(i-1,j,k,nrhs)-2.0_rt*vold(i ,j,k,nrhs)+vold(i+1,j,k,nrhs);

if (i == dlo.x and (bcr_x.lo(0) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
if (i == dlo.x and !is_periodic_in_x) {
vxx_im1 = vxx_i;
} else if (i == dhi.x+1 and (bcr_x.hi(0) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
} else if (i == dhi.x+1 and !is_periodic_in_x) {
vxx_i = vxx_im1;
}

Expand All @@ -291,11 +295,11 @@ REMORA::rhs_uv_3d (const Box& xbx, const Box& ybx,
Real Hvee_j = Hvom(i,j-1,k)-2.0_rt*Hvom(i,j ,k)+Hvom(i,j+1,k);
Real Hvee_jp1 = Hvom(i,j ,k)-2.0_rt*Hvom(i,j+1,k)+Hvom(i,j+2,k);

if (j == dlo.y and (bcr_y.lo(1) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
if (j == dlo.y and !is_periodic_in_y) {
vee_j = vee_jp1;
Hvee_j = Hvee_jp1;
}
else if (j == dhi.y and (bcr_y.hi(1) == REMORABCType::ext_dir or ic_bc_type==IC_BC_Type::Real)) {
else if (j == dhi.y and !is_periodic_in_y) {
vee_jp1 = vee_j;
Hvee_jp1 = Hvee_j;
}
Expand Down
Loading

0 comments on commit 2d2e8a4

Please sign in to comment.