Skip to content

Commit

Permalink
Merge branch 'development' into ThetavFix
Browse files Browse the repository at this point in the history
  • Loading branch information
AMLattanzi authored Dec 9, 2024
2 parents 06ecd75 + add69a0 commit f0f0af3
Show file tree
Hide file tree
Showing 2 changed files with 29 additions and 29 deletions.
31 changes: 13 additions & 18 deletions Source/IO/ERF_Write1DProfiles.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -356,33 +356,28 @@ ERF::derive_diag_profiles(Real /*time*/,
const Array4<Real>& v_cc_arr = v_cc.array(mfi);
const Array4<Real>& w_cc_arr = w_cc.array(mfi);
const Array4<Real>& p0_arr = p_hse.array(mfi);
const Array4<Real>& qv_arr = qmoist[0][0]->array(mfi); // TODO: Is this written only on lev 0?

int rhoqr_comp = solverChoice.RhoQr_comp;

ParallelFor(bx, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept
{
Real p = getPgivenRTh(cons_arr(i, j, k, RhoTheta_comp), qv_arr(i,j,k));
Real qv = cons_arr(i,j,k,RhoQ1_comp) / cons_arr(i,j,k,Rho_comp);
Real qc = cons_arr(i,j,k,RhoQ2_comp) / cons_arr(i,j,k,Rho_comp);
Real qr = (rhoqr_comp > -1) ? cons_arr(i,j,k,rhoqr_comp) / cons_arr(i,j,k,Rho_comp) :
Real(0.0);
Real p = getPgivenRTh(cons_arr(i, j, k, RhoTheta_comp), qv);

p -= p0_arr(i,j,k);
fab_arr(i, j, k,18) = p; // p
fab_arr(i, j, k,19) = p * u_cc_arr(i,j,k); // p*u
fab_arr(i, j, k,20) = p * v_cc_arr(i,j,k); // p*v
fab_arr(i, j, k,21) = p * w_cc_arr(i,j,k); // p*w
fab_arr(i, j, k,22) = cons_arr(i,j,k,RhoQ1_comp) / cons_arr(i,j,k,Rho_comp); // qv
fab_arr(i, j, k,23) = cons_arr(i,j,k,RhoQ2_comp) / cons_arr(i,j,k,Rho_comp); // qc
if (rhoqr_comp > -1) {
fab_arr(i, j, k,24) = cons_arr(i,j,k,rhoqr_comp) / cons_arr(i,j,k,Rho_comp); // qr
} else {
fab_arr(i, j, k,24) = Real(0.0);
}
fab_arr(i, j, k,25) = w_cc_arr(i,j,k) * cons_arr(i,j,k,RhoQ1_comp) / cons_arr(i,j,k,Rho_comp); // w*qv
fab_arr(i, j, k,26) = w_cc_arr(i,j,k) * cons_arr(i,j,k,RhoQ2_comp) / cons_arr(i,j,k,Rho_comp); // w*qc
if (rhoqr_comp > -1) {
fab_arr(i, j, k,27) = w_cc_arr(i,j,k) * cons_arr(i,j,k,rhoqr_comp) / cons_arr(i,j,k,Rho_comp); // w*qr
} else {
fab_arr(i, j, k,27) = Real(0.0);
}
fab_arr(i, j, k,22) = qv; // qv
fab_arr(i, j, k,23) = qc; // qc
fab_arr(i, j, k,24) = qr; // qr
fab_arr(i, j, k,25) = w_cc_arr(i,j,k) * qv; // w*qv
fab_arr(i, j, k,26) = w_cc_arr(i,j,k) * qc; // w*qc
fab_arr(i, j, k,27) = w_cc_arr(i,j,k) * qr; // w*qr
if (n_qstate > 3) {
fab_arr(i, j, k,28) = cons_arr(i,j,k,RhoQ3_comp) / cons_arr(i,j,k,Rho_comp); // qi
fab_arr(i, j, k,29) = cons_arr(i,j,k,RhoQ5_comp) / cons_arr(i,j,k,Rho_comp); // qs
Expand All @@ -392,9 +387,9 @@ ERF::derive_diag_profiles(Real /*time*/,
fab_arr(i, j, k,29) = 0.0; // qs
fab_arr(i, j, k,30) = 0.0; // qg
}
Real ql = fab_arr(i, j, k,22) + fab_arr(i, j, k,23);
Real ql = qc + qr;
Real theta = cons_arr(i,j,k,RhoTheta_comp) / cons_arr(i,j,k,Rho_comp);
Real thv = theta * (1 + 0.61*qv_arr(i,j,k) - ql);
Real thv = theta * (1 + 0.61*qv - ql);
fab_arr(i, j, k,31) = w_cc_arr(i,j,k) * thv; // w*thv
});
} // mfi
Expand Down
27 changes: 16 additions & 11 deletions Source/IO/ERF_Write1DProfiles_stag.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -459,21 +459,24 @@ ERF::derive_diag_profiles_stag (Real /*time*/,
const Array4<Real>& v_cc_arr = v_cc.array(mfi);
const Array4<Real>& w_fc_arr = w_fc.array(mfi);
const Array4<Real>& p0_arr = p_hse.array(mfi);
const Array4<Real>& qv_arr = qmoist[0][0]->array(mfi); // TODO: Is this written only on lev 0?

int rhoqr_comp = solverChoice.RhoQr_comp;

ParallelFor(bx, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept
{
Real p = getPgivenRTh(cons_arr(i, j, k, RhoTheta_comp), qv_arr(i,j,k));
Real qv = cons_arr(i,j,k,RhoQ1_comp) / cons_arr(i,j,k,Rho_comp);
Real qc = cons_arr(i,j,k,RhoQ2_comp) / cons_arr(i,j,k,Rho_comp);
Real qr = (rhoqr_comp > -1) ? cons_arr(i,j,k,rhoqr_comp) / cons_arr(i,j,k,Rho_comp) :
Real(0.0);
Real p = getPgivenRTh(cons_arr(i, j, k, RhoTheta_comp), qv);

p -= p0_arr(i,j,k);
fab_arr(i, j, k,13) = p; // p
fab_arr(i, j, k,14) = p * u_cc_arr(i,j,k); // p*u
fab_arr(i, j, k,15) = p * v_cc_arr(i,j,k); // p*v
fab_arr(i, j, k,16) = cons_arr(i,j,k,RhoQ1_comp) / cons_arr(i,j,k,Rho_comp); // qv
fab_arr(i, j, k,17) = cons_arr(i,j,k,RhoQ2_comp) / cons_arr(i,j,k,Rho_comp); // qc
fab_arr(i, j, k,18) = cons_arr(i,j,k,rhoqr_comp) / cons_arr(i,j,k,Rho_comp); // qr
fab_arr(i, j, k,16) = qv; // qv
fab_arr(i, j, k,17) = qc; // qc
fab_arr(i, j, k,18) = qr; // qr
if (n_qstate > 3) { // SAM model
fab_arr(i, j, k,19) = cons_arr(i,j,k,RhoQ3_comp) / cons_arr(i,j,k,Rho_comp); // qi
fab_arr(i, j, k,20) = cons_arr(i,j,k,RhoQ5_comp) / cons_arr(i,j,k,Rho_comp); // qs
Expand All @@ -488,20 +491,22 @@ ERF::derive_diag_profiles_stag (Real /*time*/,
const Box& zbx = mfi.tilebox(IntVect(0,0,1));
ParallelFor(zbx, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept
{
Real p0 = getPgivenRTh(cons_arr(i, j, k , RhoTheta_comp), qv_arr(i,j,k )) - p0_arr(i,j,k );
Real p1 = getPgivenRTh(cons_arr(i, j, k-1, RhoTheta_comp), qv_arr(i,j,k-1)) - p0_arr(i,j,k-1);
Real pface = 0.5 * (p0 + p1);

Real qv0 = cons_arr(i,j,k ,RhoQ1_comp) / cons_arr(i,j,k ,Rho_comp);
Real qv1 = cons_arr(i,j,k-1,RhoQ1_comp) / cons_arr(i,j,k-1,Rho_comp);
Real qc0 = cons_arr(i,j,k ,RhoQ2_comp) / cons_arr(i,j,k ,Rho_comp);
Real qc1 = cons_arr(i,j,k-1,RhoQ2_comp) / cons_arr(i,j,k-1,Rho_comp);
Real qr0 = cons_arr(i,j,k ,RhoQ3_comp) / cons_arr(i,j,k ,Rho_comp);
Real qr1 = cons_arr(i,j,k-1,RhoQ3_comp) / cons_arr(i,j,k-1,Rho_comp);
Real qr0 = (rhoqr_comp > -1) ? cons_arr(i,j,k ,RhoQ3_comp) / cons_arr(i,j,k ,Rho_comp) :
Real(0.0);
Real qr1 = (rhoqr_comp > -1) ? cons_arr(i,j,k-1,RhoQ3_comp) / cons_arr(i,j,k-1,Rho_comp) :
Real(0.0);
Real qvface = 0.5 * (qv0 + qv1);
Real qcface = 0.5 * (qc0 + qc1);
Real qrface = 0.5 * (qr0 + qr1);

Real p0 = getPgivenRTh(cons_arr(i, j, k , RhoTheta_comp), qv0) - p0_arr(i,j,k );
Real p1 = getPgivenRTh(cons_arr(i, j, k-1, RhoTheta_comp), qv1) - p0_arr(i,j,k-1);
Real pface = 0.5 * (p0 + p1);

Real theta0 = cons_arr(i,j,k ,RhoTheta_comp) / cons_arr(i,j,k ,Rho_comp);
Real theta1 = cons_arr(i,j,k-1,RhoTheta_comp) / cons_arr(i,j,k-1,Rho_comp);
Real thface = 0.5*(theta0 + theta1);
Expand Down

0 comments on commit f0f0af3

Please sign in to comment.