From e7472e2642ce1aa479368ee9b457c104dfe1e8ce Mon Sep 17 00:00:00 2001 From: japham0 Date: Fri, 21 Jun 2024 14:57:34 -0700 Subject: [PATCH 1/9] changes to read land/sea mask --- Source/ERF_read_waves.cpp | 55 ++++++++++++++++----------------------- 1 file changed, 23 insertions(+), 32 deletions(-) diff --git a/Source/ERF_read_waves.cpp b/Source/ERF_read_waves.cpp index 9af9410b4..524917c80 100644 --- a/Source/ERF_read_waves.cpp +++ b/Source/ERF_read_waves.cpp @@ -12,12 +12,10 @@ ERF::read_waves (int lev) { for ( MFIter mfi(*Hwave_onegrid[lev],false); mfi.isValid(); ++mfi) { - //auto my_H_ptr = Hwave[lev]->array(mfi); - //auto my_L_ptr = Lwave[lev]->array(mfi); + const auto & bx = mfi.validbox(); - amrex::Print() << " HERE " << bx << std::endl; - // How to declare my_H_ptr directly? + amrex::Print() << " HERE " << bx << std::endl; amrex::Array4 my_H_arr = Hwave_onegrid[lev]->array(mfi); amrex::Array4 my_L_arr = Lwave_onegrid[lev]->array(mfi); @@ -35,17 +33,17 @@ ERF::read_waves (int lev) } - int nx=2147483647; // sanity check + int nx=2147483647; int ny=2147483647; // sanity check - //JUST RECV + //JUST RECEIVED if (amrex::MPMD::MyProc() == this_root) { - if (rank_offset == 0) // the first program + if (rank_offset == 0) // First program { MPI_Recv(&nx, 1, MPI_INT, other_root, 1, MPI_COMM_WORLD, MPI_STATUS_IGNORE); MPI_Recv(&ny, 1, MPI_INT, other_root, 7, MPI_COMM_WORLD, MPI_STATUS_IGNORE); } - else // the second program + else // Second program { MPI_Recv(&nx, 1, MPI_INT, other_root, 0, MPI_COMM_WORLD, MPI_STATUS_IGNORE); MPI_Recv(&ny, 1, MPI_INT, other_root, 6, MPI_COMM_WORLD, MPI_STATUS_IGNORE); @@ -54,17 +52,9 @@ ERF::read_waves (int lev) ParallelDescriptor::Bcast(&nx, 1); ParallelDescriptor::Bcast(&ny, 1); } + if((nx)*(ny) > 0) { int nsealm = (nx)*ny; - Print()<ParallelCopy(*Hwave_onegrid[lev]); @@ -99,9 +80,19 @@ ERF::read_waves (int lev) Lwave[lev]->ParallelCopy(*Lwave_onegrid[lev]); Lwave[lev]->FillBoundary(geom[lev].periodicity()); amrex::Print() << "HWAVE BOX " << (*Hwave[lev])[0].box() << std::endl; - - print_state(*Hwave[lev],IntVect(103,-3,0),0,IntVect(3,3,0)); - print_state(*Hwave[lev],IntVect(103,88,0),0,IntVect(3,3,0)); + for (MFIter mfi(*Hwave[lev],TilingIfNotGPU()); mfi.isValid(); ++mfi) { + Box bx = mfi.tilebox(); + const Array4& Hwave_arr = Hwave[lev]->const_array(mfi); + const Array4& Lmask_arr = lmask_lev[lev][0]->array(mfi); + ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k){ + if (Hwave_arr(i,j,k)<0) { + Lmask_arr(i,j,k) = 1; + } else { + Lmask_arr(i,j,k) = 0; + } + }); + } + } #endif From a5634af674d2f6e847c5eb79dfe0c7c90ec7321c Mon Sep 17 00:00:00 2001 From: japham0 Date: Fri, 21 Jun 2024 15:10:52 -0700 Subject: [PATCH 2/9] changes to read land/sea mask --- Source/ERF_read_waves.cpp | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/Source/ERF_read_waves.cpp b/Source/ERF_read_waves.cpp index 524917c80..ec9032992 100644 --- a/Source/ERF_read_waves.cpp +++ b/Source/ERF_read_waves.cpp @@ -12,7 +12,7 @@ ERF::read_waves (int lev) { for ( MFIter mfi(*Hwave_onegrid[lev],false); mfi.isValid(); ++mfi) { - + const auto & bx = mfi.validbox(); amrex::Print() << " HERE " << bx << std::endl; @@ -33,7 +33,7 @@ ERF::read_waves (int lev) } - int nx=2147483647; + int nx=2147483647; int ny=2147483647; // sanity check //JUST RECEIVED @@ -69,7 +69,7 @@ ERF::read_waves (int lev) } } - } + } } @@ -92,7 +92,7 @@ ERF::read_waves (int lev) } }); } - + } #endif From b77d07dbe2e7e740a91df8f99578487913201588 Mon Sep 17 00:00:00 2001 From: japham0 Date: Wed, 26 Jun 2024 14:06:09 -0700 Subject: [PATCH 3/9] send_waves function to send velocities to WW3 --- Source/ERF.H | 1 + Source/ERF.cpp | 1 + Source/ERF_read_waves.cpp | 81 +++++++++++++++++++++++++ Source/TimeIntegration/ERF_TimeStep.cpp | 1 + 4 files changed, 84 insertions(+) diff --git a/Source/ERF.H b/Source/ERF.H index 4cb59e790..aacef2395 100644 --- a/Source/ERF.H +++ b/Source/ERF.H @@ -240,6 +240,7 @@ public: #ifdef ERF_USE_WW3_COUPLING void read_waves (int lev); + void send_waves (int lev); #endif // Interface for advancing the data at one level by one "slow" timestep diff --git a/Source/ERF.cpp b/Source/ERF.cpp index 26faab742..b5ff512bb 100644 --- a/Source/ERF.cpp +++ b/Source/ERF.cpp @@ -909,6 +909,7 @@ ERF::InitData () #ifdef ERF_USE_WW3_COUPLING int lev = 0; read_waves(lev); + send_waves(lev); #endif // Configure ABLMost params if used MostWall boundary condition diff --git a/Source/ERF_read_waves.cpp b/Source/ERF_read_waves.cpp index ec9032992..a00fc32f6 100644 --- a/Source/ERF_read_waves.cpp +++ b/Source/ERF_read_waves.cpp @@ -4,6 +4,8 @@ #include #include #include +#include +#include using namespace amrex; @@ -68,6 +70,8 @@ ERF::read_waves (int lev) MPI_Recv(my_L_ptr, nsealm, MPI_DOUBLE, other_root, 4, MPI_COMM_WORLD, MPI_STATUS_IGNORE); } } + amrex::AllPrintToFile("output_HS_cpp.txt")<ParallelCopy(*Lwave_onegrid[lev]); Lwave[lev]->FillBoundary(geom[lev].periodicity()); amrex::Print() << "HWAVE BOX " << (*Hwave[lev])[0].box() << std::endl; + for (MFIter mfi(*Hwave[lev],TilingIfNotGPU()); mfi.isValid(); ++mfi) { Box bx = mfi.tilebox(); const Array4& Hwave_arr = Hwave[lev]->const_array(mfi); @@ -93,6 +98,82 @@ ERF::read_waves (int lev) }); } +} + +void +ERF::send_waves (int lev) +{ + int ncomp = 1; // number components + auto& lev_new = vars_new[lev]; + + // Access xvel, yvel from ABL + MultiFab xvel_data(lev_new[Vars::xvel].boxArray(), lev_new[Vars::xvel].DistributionMap(), 1, lev_new[Vars::xvel].nGrowVect()); + MultiFab yvel_data(lev_new[Vars::yvel].boxArray(), lev_new[Vars::yvel].DistributionMap(), 1, lev_new[Vars::yvel].nGrowVect()); + + // Make local copy of xvel, yvel + MultiFab::Copy (xvel_data, lev_new[Vars::xvel], 0, 0, 1, lev_new[Vars::xvel].nGrowVect()); + MultiFab::Copy (yvel_data, lev_new[Vars::yvel], 0, 0, 1, lev_new[Vars::yvel].nGrowVect()); + + + MultiFab x_avg(lev_new[Vars::cons].boxArray(), lev_new[Vars::cons].DistributionMap(), + ncomp, lev_new[Vars::cons].nGrow()); + MultiFab y_avg(lev_new[Vars::cons].boxArray(), lev_new[Vars::cons].DistributionMap(), + ncomp, lev_new[Vars::cons].nGrow()); + + MultiFab u_mag(lev_new[Vars::cons].boxArray(), lev_new[Vars::cons].DistributionMap(), + ncomp, lev_new[Vars::cons].nGrow()); + + MultiFab u_dir(lev_new[Vars::cons].boxArray(), lev_new[Vars::cons].DistributionMap(), + ncomp, lev_new[Vars::cons].nGrow()); + x_avg.setVal(0.); + y_avg.setVal(0.); + u_mag.setVal(0.); + u_dir.setVal(0.); + + for (MFIter mfi(x_avg,TilingIfNotGPU()); mfi.isValid(); ++mfi) { + + Box bx = mfi.tilebox(); + const Array4& u_vel = x_avg.array(mfi); + const Array4& velx_arr = xvel_data.array(mfi); + + ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k){ + u_vel(i,j,k) = 0.5 *( velx_arr(i,j,k) + velx_arr(i+1,j,k) ); + + amrex::AllPrintToFile("uvel.txt") << amrex::IntVect(i,j,k) << " [" <& v_vel = y_avg.array(mfi); + const Array4& vely_arr = yvel_data.array(mfi); + + ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k){ + + v_vel(i,j,k) = 0.5 *( vely_arr(i,j,k) + vely_arr(i,j+1,k) ); + + amrex::AllPrintToFile("vvel.txt") << amrex::IntVect(i,j,k) << " ["<& magnitude = u_mag.array(mfi); + const Array4& u = x_avg.array(mfi); + const Array4& v = y_avg.array(mfi); + const Array4& theta = u_dir.array(mfi); + + ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k){ + + magnitude(i,j,k) = std::sqrt( pow(u(i,j,k), 2) + pow(v(i,j,k), 2) ); + theta(i,j,k) = atan ( v(i,j,k) / u(i,j,k) ); + + amrex::AllPrintToFile("mag_theta.txt") << amrex::IntVect(i,j,k) << " Magnitude: " << magnitude(i,j,k) << " Theta: " << theta(i,j,k) < Date: Mon, 1 Jul 2024 14:57:48 -0700 Subject: [PATCH 4/9] Updates to computing theta and documentation --- Docs/sphinx_doc/index.rst | 1 + Source/ERF_read_waves.cpp | 28 +++++++++++++++++++++------- 2 files changed, 22 insertions(+), 7 deletions(-) diff --git a/Docs/sphinx_doc/index.rst b/Docs/sphinx_doc/index.rst index 178d0b7fe..6931338ae 100644 --- a/Docs/sphinx_doc/index.rst +++ b/Docs/sphinx_doc/index.rst @@ -79,6 +79,7 @@ In addition to this documentation, there is API documentation for ERF generated CouplingToAMRWind.rst + CouplingToWW3.rst .. toctree:: :caption: ERF vs WRF :maxdepth: 1 diff --git a/Source/ERF_read_waves.cpp b/Source/ERF_read_waves.cpp index a00fc32f6..66c4854ea 100644 --- a/Source/ERF_read_waves.cpp +++ b/Source/ERF_read_waves.cpp @@ -105,16 +105,17 @@ ERF::send_waves (int lev) { int ncomp = 1; // number components auto& lev_new = vars_new[lev]; + const double PI = 3.1415926535897932384626433832795028841971693993751058209; // Access xvel, yvel from ABL MultiFab xvel_data(lev_new[Vars::xvel].boxArray(), lev_new[Vars::xvel].DistributionMap(), 1, lev_new[Vars::xvel].nGrowVect()); MultiFab yvel_data(lev_new[Vars::yvel].boxArray(), lev_new[Vars::yvel].DistributionMap(), 1, lev_new[Vars::yvel].nGrowVect()); - + // Make local copy of xvel, yvel MultiFab::Copy (xvel_data, lev_new[Vars::xvel], 0, 0, 1, lev_new[Vars::xvel].nGrowVect()); MultiFab::Copy (yvel_data, lev_new[Vars::yvel], 0, 0, 1, lev_new[Vars::yvel].nGrowVect()); - + MultiFab x_avg(lev_new[Vars::cons].boxArray(), lev_new[Vars::cons].DistributionMap(), ncomp, lev_new[Vars::cons].nGrow()); MultiFab y_avg(lev_new[Vars::cons].boxArray(), lev_new[Vars::cons].DistributionMap(), @@ -133,7 +134,7 @@ ERF::send_waves (int lev) for (MFIter mfi(x_avg,TilingIfNotGPU()); mfi.isValid(); ++mfi) { Box bx = mfi.tilebox(); - const Array4& u_vel = x_avg.array(mfi); + const Array4& u_vel = x_avg.array(mfi); const Array4& velx_arr = xvel_data.array(mfi); ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k){ @@ -161,19 +162,32 @@ ERF::send_waves (int lev) Box bx = mfi.tilebox(); const Array4& magnitude = u_mag.array(mfi); - const Array4& u = x_avg.array(mfi); - const Array4& v = y_avg.array(mfi); + const Array4& u = x_avg.array(mfi); + const Array4& v = y_avg.array(mfi); const Array4& theta = u_dir.array(mfi); ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k){ magnitude(i,j,k) = std::sqrt( pow(u(i,j,k), 2) + pow(v(i,j,k), 2) ); - theta(i,j,k) = atan ( v(i,j,k) / u(i,j,k) ); + + if ( u(i,j,k) < 0 && v(i,j,k) > 0 || u(i,j,k) < 0 && v(i,j,k) < 0 ) { + + theta(i,j,k) = PI + ( atan( v(i,j,k) / u(i,j,k) ) ); + + } else { + + theta(i,j,k) = atan ( v(i,j,k) / u(i,j,k) ); + } + amrex::AllPrintToFile("mag_theta.txt") << amrex::IntVect(i,j,k) << " Magnitude: " << magnitude(i,j,k) << " Theta: " << theta(i,j,k) < Date: Mon, 1 Jul 2024 15:12:13 -0700 Subject: [PATCH 5/9] try to fix .rst file --- Docs/sphinx_doc/CouplingToWW3.rst | 34 +++++++++++++++++++++++++++++++ Docs/sphinx_doc/index.rst | 5 +++++ 2 files changed, 39 insertions(+) create mode 100644 Docs/sphinx_doc/CouplingToWW3.rst diff --git a/Docs/sphinx_doc/CouplingToWW3.rst b/Docs/sphinx_doc/CouplingToWW3.rst new file mode 100644 index 000000000..7c6ed7438 --- /dev/null +++ b/Docs/sphinx_doc/CouplingToWW3.rst @@ -0,0 +1,34 @@ + + .. role:: cpp(code) + :language: c++ + + .. _CouplingToWW3: + +Coupling To WW3 +=============== + +Coupling with WaveWatch III is currently a work in progress. +Currently, we have a one-way coupling between ERF and WaveWatch III (WW3), where WW3 sends ERF Hwave (significant wave height) and Lwave (mean wavelength) over a grid. The + +One-way coupling WW3 to ERF +--------------------------- + +Values are used to compute the surface roughness z0 through a fixed-point iteration: + +.. math:: + + z0 = 1200.0 Hwave \left(\frac{Hwave}{Lwave}\right)^{4.5} + \frac{0.11 \mu}{u_*} + +To run the coupled model: + +.. code-block:: bash + + git clone --recursive git@github.com:erf-model/ERF + cd ERF/Exec/ABL + make -j4 USE_WW3_COUPLING=TRUE + cd ../../Submodules/WW3 + ./model/bin/w3_setup model -c gnu -s Ifremer1 + cd regtests + ./bin/run_cmake_test -C MPMD -n 2 -p mpirun -f -s PR1_MPI ../model ww3_tp2.2 + + diff --git a/Docs/sphinx_doc/index.rst b/Docs/sphinx_doc/index.rst index 6931338ae..7f73dec2a 100644 --- a/Docs/sphinx_doc/index.rst +++ b/Docs/sphinx_doc/index.rst @@ -79,7 +79,12 @@ In addition to this documentation, there is API documentation for ERF generated CouplingToAMRWind.rst +.. toctree:: + :caption: ERF vs WRF + :maxdepth: 1 + CouplingToWW3.rst + .. toctree:: :caption: ERF vs WRF :maxdepth: 1 From 35e87e42621ee0a7ae9e89bca52bc9e00e4c002a Mon Sep 17 00:00:00 2001 From: japham0 Date: Tue, 2 Jul 2024 10:26:04 -0700 Subject: [PATCH 6/9] updates to ERF_send_waves and docs --- Docs/sphinx_doc/CouplingToWW3.rst | 20 ++++++++++++++ Source/ERF_read_waves.cpp | 46 +++++++++++++++++++++++++++++-- 2 files changed, 63 insertions(+), 3 deletions(-) diff --git a/Docs/sphinx_doc/CouplingToWW3.rst b/Docs/sphinx_doc/CouplingToWW3.rst index 7c6ed7438..2b98a59dd 100644 --- a/Docs/sphinx_doc/CouplingToWW3.rst +++ b/Docs/sphinx_doc/CouplingToWW3.rst @@ -31,4 +31,24 @@ To run the coupled model: cd regtests ./bin/run_cmake_test -C MPMD -n 2 -p mpirun -f -s PR1_MPI ../model ww3_tp2.2 +Modifications to the problem size and geometry, as well as other parameters can be done in the `inputs_mpmd` file. The `plt` files as well as relevant outputs can be viewed in the `regtests/ww3_tp2.2/work` directory. + +Two-way coupling: +----------------- + +Disclaimer: Two-way coupling is currently a work in progress. Two-way coupling involves sending the wind velocity and direction to WW3. We convert the x and y velocities from ERF to a wind speed (U_{wind}) and wind direction (\theta). + +.. math:: + + \text{Wind speed at reference height} \quad U_{wind} = \sqrt{u^2 + v^2} + +.. math:: + + \text{Wind direction at reference height} \quad \theta = \arctan{\frac{v}{u}} + +Both U_{wind} and \theta are then used in the wind source term S_{in} in the ST6 subroutine in WW3. To run the model with two-way coupling: + +.. code-block:: bash + + ./bin/run_cmake_test -C MPMD -n 2 -p mpirun -f -s PR1_ST6 ../model ww3_ts2 diff --git a/Source/ERF_read_waves.cpp b/Source/ERF_read_waves.cpp index 66c4854ea..2d45d5a7d 100644 --- a/Source/ERF_read_waves.cpp +++ b/Source/ERF_read_waves.cpp @@ -170,6 +170,13 @@ ERF::send_waves (int lev) magnitude(i,j,k) = std::sqrt( pow(u(i,j,k), 2) + pow(v(i,j,k), 2) ); + double u_val = u(i, j, k); + + if ( u_val == 0 ) { + u_val = std::max( u_val, 1e-15 ); // Ensure u_val is non-zero + } + + if ( u(i,j,k) < 0 && v(i,j,k) > 0 || u(i,j,k) < 0 && v(i,j,k) < 0 ) { theta(i,j,k) = PI + ( atan( v(i,j,k) / u(i,j,k) ) ); @@ -184,10 +191,43 @@ ERF::send_waves (int lev) }); + // MPI_Send to WW3 + // Calculate the number of elements in the current box + int n_elements = bx.numPts(); - //Real* magnitude = magnitude.dataPtr(); - //Real* theta = theta.dataPtr(); - } + // Get the data pointers for the arrays + Real* magnitude_ptr = &magnitude(bx.smallEnd()); + Real* theta_ptr = &theta(bx.smallEnd()); + + // Initialize other_root as needed + int other_root = 0; // Example initialization, replace with appropriate logic + + // Send magnitude array +// MPI_Send(magnitude_ptr, n_elements, MPI_DOUBLE, other_root, 0, MPI_COMM_WORLD); + // Send theta array +// MPI_Send(theta_ptr, n_elements, MPI_DOUBLE, other_root, 1, MPI_COMM_WORLD); +amrex::AllPrintToFile("debug_send.txt") << n_elements << std::endl; + +} +/* + for (MFIter mfi(u_mag); mfi.isValid(); ++mfi) { + + const Array4& magnitude = u_mag.array(mfi); + const Array4& theta = u_dir.array(mfi); + + Real* magnitude_ptr = magnitude.dataPtr(); + Real* theta_ptr = theta.dataPtr(); + + // Get number of elements in arrays + int n_elements = mfi.validbox().numPts(); + int this_root = 0; + int other_root = 1; + + // Send magnitude and theta + MPI_Send(magnitude_ptr. n_elements, MPI_DOUBLE, this_root, 0, MPI_COMM_WORLD) + MPI_Send(theta_ptr, n_elements, MPI_DOUBLE, other_root, 1, MPI_COMM_WORLD); + } +*/ } #endif From e641672af460d4b9a220858914f344c411825046 Mon Sep 17 00:00:00 2001 From: japham0 Date: Tue, 2 Jul 2024 11:30:31 -0700 Subject: [PATCH 7/9] fix to doc --- Docs/sphinx_doc/CouplingToWW3.rst | 20 +++++++------------- 1 file changed, 7 insertions(+), 13 deletions(-) diff --git a/Docs/sphinx_doc/CouplingToWW3.rst b/Docs/sphinx_doc/CouplingToWW3.rst index 2b98a59dd..5f8e9b5d9 100644 --- a/Docs/sphinx_doc/CouplingToWW3.rst +++ b/Docs/sphinx_doc/CouplingToWW3.rst @@ -8,16 +8,15 @@ Coupling To WW3 =============== Coupling with WaveWatch III is currently a work in progress. -Currently, we have a one-way coupling between ERF and WaveWatch III (WW3), where WW3 sends ERF Hwave (significant wave height) and Lwave (mean wavelength) over a grid. The +Currently, we have a one-way coupling between ERF and WaveWatch III (WW3), where WW3 sends ERF Hwave (significant wave height) and Lwave (mean wavelength) over a grid. One-way coupling WW3 to ERF --------------------------- -Values are used to compute the surface roughness z0 through a fixed-point iteration: +The values are used to compute the surface roughness :math:`\overline{z_{0}}` through a fixed-point iteration: .. math:: - - z0 = 1200.0 Hwave \left(\frac{Hwave}{Lwave}\right)^{4.5} + \frac{0.11 \mu}{u_*} + \overline{z_{0}} = 1200.0 Hwave (\frac{Hwave}{Lwave})^{4.5} + \frac{0.11 \mu}{u_*} To run the coupled model: @@ -36,19 +35,14 @@ Modifications to the problem size and geometry, as well as other parameters can Two-way coupling: ----------------- -Disclaimer: Two-way coupling is currently a work in progress. Two-way coupling involves sending the wind velocity and direction to WW3. We convert the x and y velocities from ERF to a wind speed (U_{wind}) and wind direction (\theta). +Disclaimer: Two-way coupling is currently a work in progress. Two-way coupling involves sending the wind velocity and direction to WW3. We convert the x and y velocities from ERF to a wind speed (:math:`\overline{U_{wind}}`) and wind direction (:math:`\theta`) from the reference height. .. math:: - \text{Wind speed at reference height} \quad U_{wind} = \sqrt{u^2 + v^2} + \overline{U_{wind}} = \sqrt{\overline{u^{2}} + \overline{v^2}} .. math:: - \text{Wind direction at reference height} \quad \theta = \arctan{\frac{v}{u}} - -Both U_{wind} and \theta are then used in the wind source term S_{in} in the ST6 subroutine in WW3. To run the model with two-way coupling: - -.. code-block:: bash - - ./bin/run_cmake_test -C MPMD -n 2 -p mpirun -f -s PR1_ST6 ../model ww3_ts2 + \overline{\theta} = \mathrm{arctan}{\frac{v}{u}} +Both :math:`\overline{U_{wind}}` and :math:`\theta` are then used in the wind source term :math:`S_{in}` in the ST6 subroutine in WW3 From f953b42218e188afe2833c0ecbd050d6d43c1696 Mon Sep 17 00:00:00 2001 From: japham0 Date: Tue, 2 Jul 2024 11:41:12 -0700 Subject: [PATCH 8/9] fixes to doc --- Docs/sphinx_doc/CouplingToWW3.rst | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/Docs/sphinx_doc/CouplingToWW3.rst b/Docs/sphinx_doc/CouplingToWW3.rst index 5f8e9b5d9..cd648b260 100644 --- a/Docs/sphinx_doc/CouplingToWW3.rst +++ b/Docs/sphinx_doc/CouplingToWW3.rst @@ -16,7 +16,7 @@ One-way coupling WW3 to ERF The values are used to compute the surface roughness :math:`\overline{z_{0}}` through a fixed-point iteration: .. math:: - \overline{z_{0}} = 1200.0 Hwave (\frac{Hwave}{Lwave})^{4.5} + \frac{0.11 \mu}{u_*} + z_{0} = 1200.0 Hwave (\frac{Hwave}{Lwave})^{4.5} + \frac{0.11 \mu}{u_*} To run the coupled model: @@ -35,14 +35,14 @@ Modifications to the problem size and geometry, as well as other parameters can Two-way coupling: ----------------- -Disclaimer: Two-way coupling is currently a work in progress. Two-way coupling involves sending the wind velocity and direction to WW3. We convert the x and y velocities from ERF to a wind speed (:math:`\overline{U_{wind}}`) and wind direction (:math:`\theta`) from the reference height. +Disclaimer: Two-way coupling is currently a work in progress. Two-way coupling involves sending the wind velocity and direction to WW3. We convert the x and y velocities from ERF to a wind speed (:math:`U_{wind}`) and wind direction (:math:`\theta`) from the reference height. .. math:: - \overline{U_{wind}} = \sqrt{\overline{u^{2}} + \overline{v^2}} + U_{wind} = \sqrt{u^{2} + v^2} .. math:: - \overline{\theta} = \mathrm{arctan}{\frac{v}{u}} + \theta = \mathrm{arctan}{\frac{v}{u}} -Both :math:`\overline{U_{wind}}` and :math:`\theta` are then used in the wind source term :math:`S_{in}` in the ST6 subroutine in WW3 +Both :math:`U_{wind}` and :math:`\theta` are then used in the wind source term :math:`S_{in}` in the ST6 subroutine in WW3 From 77971236b5e51124b660b06821653604a183d769 Mon Sep 17 00:00:00 2001 From: japham0 Date: Tue, 2 Jul 2024 11:46:52 -0700 Subject: [PATCH 9/9] fixed whitespaces --- Docs/sphinx_doc/CouplingToWW3.rst | 2 +- Exec/ABL/inputs_mpmd | 4 ++-- Source/ERF_read_waves.cpp | 8 ++++---- 3 files changed, 7 insertions(+), 7 deletions(-) diff --git a/Docs/sphinx_doc/CouplingToWW3.rst b/Docs/sphinx_doc/CouplingToWW3.rst index cd648b260..5a9934d0c 100644 --- a/Docs/sphinx_doc/CouplingToWW3.rst +++ b/Docs/sphinx_doc/CouplingToWW3.rst @@ -8,7 +8,7 @@ Coupling To WW3 =============== Coupling with WaveWatch III is currently a work in progress. -Currently, we have a one-way coupling between ERF and WaveWatch III (WW3), where WW3 sends ERF Hwave (significant wave height) and Lwave (mean wavelength) over a grid. +Currently, we have a one-way coupling between ERF and WaveWatch III (WW3), where WW3 sends ERF Hwave (significant wave height) and Lwave (mean wavelength) over a grid. One-way coupling WW3 to ERF --------------------------- diff --git a/Exec/ABL/inputs_mpmd b/Exec/ABL/inputs_mpmd index d03df316a..a4ffbd1dc 100644 --- a/Exec/ABL/inputs_mpmd +++ b/Exec/ABL/inputs_mpmd @@ -6,8 +6,8 @@ amrex.fpe_trap_invalid = 0 fabarray.mfiter_tile_size = 1024 1024 1024 # PROBLEM SIZE & GEOMETRY -geometry.prob_extent = 191000 91000 1024 -amr.n_cell = 191 91 4 +geometry.prob_extent = 191 91 1024 +amr.n_cell = 191 91 4 geometry.is_periodic = 1 1 0 diff --git a/Source/ERF_read_waves.cpp b/Source/ERF_read_waves.cpp index 2d45d5a7d..22af1b575 100644 --- a/Source/ERF_read_waves.cpp +++ b/Source/ERF_read_waves.cpp @@ -171,19 +171,19 @@ ERF::send_waves (int lev) magnitude(i,j,k) = std::sqrt( pow(u(i,j,k), 2) + pow(v(i,j,k), 2) ); double u_val = u(i, j, k); - + double v_val = v(i, j, k); if ( u_val == 0 ) { u_val = std::max( u_val, 1e-15 ); // Ensure u_val is non-zero } - if ( u(i,j,k) < 0 && v(i,j,k) > 0 || u(i,j,k) < 0 && v(i,j,k) < 0 ) { + if ( u_val < 0 && v_val > 0 || u_val < 0 && v_val < 0 ) { - theta(i,j,k) = PI + ( atan( v(i,j,k) / u(i,j,k) ) ); + theta(i,j,k) = PI + ( atan( v_val / u_val ) ); } else { - theta(i,j,k) = atan ( v(i,j,k) / u(i,j,k) ); + theta(i,j,k) = atan ( v_val / u_val ); }