Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Mc optically thick sphere #6403

Open
wants to merge 1 commit into
base: develop
Choose a base branch
from

Conversation

ffoucart
Copy link
Contributor

Proposed changes

Add initial conditions for a homogeneous sphere (fluid variables only determined by whether we are inside or outside the sphere)

Upgrade instructions

None

Code review checklist

  • The code is documented and the documentation renders correctly. Run
    make doc to generate the documentation locally into BUILD_DIR/docs/html.
    Then open index.html.
  • The code follows the stylistic and code quality guidelines listed in the
    code review guide.
  • The PR lists upgrade instructions and is labeled bugfix or
    new feature if appropriate.

@ffoucart
Copy link
Contributor Author

Depends on #6355 (first five commits are from that PR)

@ffoucart ffoucart force-pushed the MCOpticallyThickSphereID branch from 7d7ae5f to 9bbe256 Compare December 26, 2024 18:23
@ffoucart
Copy link
Contributor Author

Updated so that first two commits match #6355 and third commit is the new work.
Note that in addition of creating input file options for MC, this PR fixes a bug in the initialization of the MortarData.

@ffoucart ffoucart mentioned this pull request Jan 24, 2025
3 tasks
@ffoucart ffoucart force-pushed the MCOpticallyThickSphereID branch from 9bbe256 to dad3fdf Compare January 27, 2025 17:36
@ffoucart ffoucart force-pushed the MCOpticallyThickSphereID branch from dad3fdf to 26b2fe5 Compare January 29, 2025 17:51
@ffoucart ffoucart marked this pull request as ready for review January 29, 2025 17:52
Copy link
Member

@kidder kidder left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

please make changes in a fixup commit

@@ -131,12 +137,20 @@ struct EvolutionMetavars {
};

using observe_fields =
tmpl::list<domain::Tags::Coordinates<volume_dim, Frame::Grid>,
domain::Tags::Coordinates<volume_dim, Frame::Inertial>>;
tmpl::list<evolution::dg::subcell::Tags::ObserverCoordinatesCompute<
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

you should make the tags of ObserverCoordinates conditional on use_dg_subcell like is done in EvolveValenciaDivClean

tmpl::list<::Events::Tags::ObserverMeshCompute<volume_dim>,
::Events::Tags::ObserverDetInvJacobianCompute<
Frame::ElementLogical, Frame::Inertial>>;
tmpl::list<evolution::dg::subcell::Tags::ObserverMeshCompute<volume_dim>,
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

again conditional on use_dg_subcell

@@ -25,6 +25,7 @@ void cell_light_crossing_time(
inertial_coordinates.get(1)[step[1]] - inertial_coordinates.get(1)[0],
inertial_coordinates.get(2)[step[2]] - inertial_coordinates.get(2)[0]};

get(*cell_light_crossing_time) = DataVector(n_pts,0.0);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

just need number of points, don't have to initialize to a value as it is overwritten

@@ -381,6 +397,26 @@ struct ReceiveDataForMcCommunication {
} else {
const size_t n_packets = received_data_packets.value().size();
for (size_t p = 0; p < n_packets; p++) {
// Find closest grid point to received packet at current
// time, using extents for live points only.
Packet& packet = received_data_packets.value()[p];
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think this can be simplified a lot as the computation is done in ElementLogicalCoordinates.
In each dimension there are n = extents[d] subcells of uniform size and the grid points are at cell centers. Therefore you just need to determine which cell they are in which is given by cell_index[d] = std::floor(packet.coordinate[d] + 1.0)/(2*n) and then index_of_closest_grid_point = collapsed_index(cell_index, extents). I guess you need to add a special case if the packet coordinate is 1.0


#include "Evolution/Particles/MonteCarlo/MonteCarloOptions.hpp"

#include <cmath>
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

don't need


DataVector get_mask(const tnsr::I<DataVector, 3>& x,
const double radius_bound) {
const DataVector radius = sqrt((x.get(0) * x.get(0)) + (x.get(1) * x.get(1)) +
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

here too use Euclidean magnitude function

@@ -65,4 +66,22 @@ struct InteractionRatesTable : db::SimpleTag {
}
};

template <size_t NeutrinoSpecies>
struct MonteCarloOptions : db::SimpleTag {
using type = std::unique_ptr<
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I don't think this needs to be a unique_ptr, but just the class

const std::vector<Packet>& packets, const Scalar<DataVector>& lapse,
const Scalar<DataVector>& sqrt_determinant_spatial_metric,
const Mesh<3>& mesh,
const Scalar<DataVector>& det_inv_jacobian_logical_to_inertial) {
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

move the function definition to the cpp file

Scalar<DataVector> cell_inertial_three_volume(lapse);
cell_inertial_coordinate_three_volume_finite_difference(
&cell_inertial_three_volume, mesh, det_jacobian_logical_to_inertial);
get(*fluid_frame_energy_density) = 0.0 * get(lapse);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

just use *fluid_frame_energy_density = make_with_value<Scalar<DataVector>>(lapse, 0.0);

using MortarData =
typename Particles::MonteCarlo::Tags::MortarDataTag<dim>::type;
MortarData mortar_data;
const Element<dim>& element = db::get<::domain::Tags::Element<dim>>(box);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do you need to do this initialization? please check

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

2 participants