Skip to content

Commit 3e87615

Browse files
committed
++
1 parent 92c0404 commit 3e87615

File tree

7 files changed

+454
-327
lines changed

7 files changed

+454
-327
lines changed

src/core/data/grid/gridlayout.hpp

Lines changed: 23 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1033,6 +1033,29 @@ namespace core
10331033

10341034

10351035

1036+
/**
1037+
* @brief BxToMoments return the indexes and associated coef to compute the linear
1038+
* interpolation necessary to project Bx onto moments.
1039+
*/
1040+
NO_DISCARD auto static constexpr BxToMoments() { return GridLayoutImpl::BxToMoments(); }
1041+
1042+
1043+
/**
1044+
* @brief ByToMoments return the indexes and associated coef to compute the linear
1045+
* interpolation necessary to project By onto moments.
1046+
*/
1047+
NO_DISCARD auto static constexpr ByToMoments() { return GridLayoutImpl::ByToMoments(); }
1048+
1049+
1050+
/**
1051+
* @brief BzToMoments return the indexes and associated coef to compute the linear
1052+
* interpolation necessary to project Bz onto moments.
1053+
*/
1054+
NO_DISCARD auto static constexpr BzToMoments() { return GridLayoutImpl::BzToMoments(); }
1055+
1056+
1057+
1058+
10361059
/**
10371060
* @brief ExToMoments return the indexes and associated coef to compute the linear
10381061
* interpolation necessary to project Ex onto moments.

src/core/data/grid/gridlayoutimplyee.hpp

Lines changed: 90 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -508,6 +508,96 @@ namespace core
508508

509509

510510

511+
NO_DISCARD auto static constexpr BxToMoments()
512+
{
513+
// Bx is primal dual dual
514+
// moments are primal primal primal
515+
// operation is thus Pdd to Ppp
516+
[[maybe_unused]] auto constexpr iShift = dualToPrimal();
517+
518+
if constexpr (dimension == 1)
519+
{
520+
constexpr WeightPoint<dimension> P1{Point<int, dimension>{0}, 1.};
521+
return std::array<WeightPoint<dimension>, 1>{P1};
522+
}
523+
if constexpr (dimension == 2)
524+
{
525+
constexpr WeightPoint<dimension> P1{Point<int, dimension>{0, 0}, 0.5};
526+
constexpr WeightPoint<dimension> P2{Point<int, dimension>{0, iShift}, 0.5};
527+
return std::array<WeightPoint<dimension>, 2>{P1, P2};
528+
}
529+
else if constexpr (dimension == 3)
530+
{
531+
constexpr WeightPoint<dimension> P1{Point<int, dimension>{0, 0, 0}, 0.5};
532+
constexpr WeightPoint<dimension> P2{Point<int, dimension>{0, iShift, iShift}, 0.5};
533+
return std::array<WeightPoint<dimension>, 2>{P1, P2};
534+
}
535+
}
536+
537+
538+
539+
540+
NO_DISCARD auto static constexpr ByToMoments()
541+
{
542+
// Bx is dual primal dual
543+
// moments are primal primal primal
544+
// operation is thus Dpd to Ppp
545+
546+
auto constexpr iShift = dualToPrimal();
547+
548+
if constexpr (dimension == 1)
549+
{
550+
constexpr WeightPoint<dimension> P1{Point<int, dimension>{0}, 0.5};
551+
constexpr WeightPoint<dimension> P2{Point<int, dimension>{iShift}, 0.5};
552+
return std::array<WeightPoint<dimension>, 2>{P1, P2};
553+
}
554+
if constexpr (dimension == 2)
555+
{
556+
constexpr WeightPoint<dimension> P1{Point<int, dimension>{0, 0}, 0.5};
557+
constexpr WeightPoint<dimension> P2{Point<int, dimension>{iShift, 0}, 0.5};
558+
return std::array<WeightPoint<dimension>, 2>{P1, P2};
559+
}
560+
else if constexpr (dimension == 3)
561+
{
562+
constexpr WeightPoint<dimension> P1{Point<int, dimension>{0, 0, 0}, 0.5};
563+
constexpr WeightPoint<dimension> P2{Point<int, dimension>{iShift, 0, iShift}, 0.5};
564+
return std::array<WeightPoint<dimension>, 2>{P1, P2};
565+
}
566+
}
567+
568+
569+
570+
571+
NO_DISCARD auto static constexpr BzToMoments()
572+
{
573+
// Bx is dual dual primal
574+
// moments are primal primal primal
575+
// operation is thus Ddp to Ppp
576+
577+
auto constexpr iShift = dualToPrimal();
578+
579+
if constexpr (dimension == 1)
580+
{
581+
constexpr WeightPoint<dimension> P1{Point<int, dimension>{0}, 0.5};
582+
constexpr WeightPoint<dimension> P2{Point<int, dimension>{iShift}, 0.5};
583+
return std::array<WeightPoint<dimension>, 2>{P1, P2};
584+
}
585+
if constexpr (dimension == 2)
586+
{
587+
constexpr WeightPoint<dimension> P1{Point<int, dimension>{0, 0}, 0.5};
588+
constexpr WeightPoint<dimension> P2{Point<int, dimension>{iShift, iShift}, 0.5};
589+
return std::array<WeightPoint<dimension>, 2>{P1, P2};
590+
}
591+
else if constexpr (dimension == 3)
592+
{
593+
constexpr WeightPoint<dimension> P1{Point<int, dimension>{0, 0, 0}, 0.5};
594+
constexpr WeightPoint<dimension> P2{Point<int, dimension>{iShift, iShift, 0}, 0.5};
595+
return std::array<WeightPoint<dimension>, 2>{P1, P2};
596+
}
597+
}
598+
599+
600+
511601

512602
NO_DISCARD auto static constexpr ExToMoments()
513603
{

src/core/utilities/algorithm.hpp

Lines changed: 17 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -77,41 +77,47 @@ auto convert_to_primal( //
7777
)
7878
{
7979
using PQ = PhysicalQuantity;
80-
constexpr std::array valid_quantities{PQ::Bx, PQ::By, PQ::Bz, PQ::Ex, PQ::Ey, PQ::Ez};
8180

82-
if (qty == PQ::Ex)
81+
if (qty == PQ::Bx)
82+
return layout.project(src, lix, layout.BxToMoments());
83+
else if (qty == PQ::By)
84+
return layout.project(src, lix, layout.ByToMoments());
85+
else if (qty == PQ::Bz)
86+
return layout.project(src, lix, layout.BzToMoments());
87+
88+
else if (qty == PQ::Ex)
8389
return layout.project(src, lix, layout.ExToMoments());
8490
else if (qty == PQ::Ey)
8591
return layout.project(src, lix, layout.EyToMoments());
8692
else if (qty == PQ::Ez)
8793
return layout.project(src, lix, layout.EzToMoments());
8894

89-
return 1e-5; // todo
95+
throw std::runtime_error("Quantity not supported for conversion to primal.");
9096
}
9197

9298
template<std::size_t dim, typename... Ts>
9399
auto& _convert_to_fortran_primal( // DOES NOT WORK ON GHOST BOX!
94-
Field<dim, Ts...>& dst, // TEMPORARY TYPE - Quantity not valid!
100+
Field<dim, Ts...>& dst, //
95101
Field<dim, Ts...> const& src, //
96102
auto const& layout //
97103
)
98104
{
99105
bool static constexpr c_ordering = false;
100-
auto const qty = src.physicalQuantity();
101-
dst.setShape(src.shape());
102106

103-
auto lb_view = core::make_array_view<c_ordering>(dst.data(), src.shape());
107+
assert(all(layout.centering(dst), [](auto const c) { return c == QtyCentering::primal; }));
108+
109+
auto lb_view = core::make_array_view<c_ordering>(dst.data(), dst.shape());
104110
auto const all_primal
105111
= all(layout.centering(src), [](auto const c) { return c == QtyCentering::primal; });
106-
auto const& amr_gbox = layout.AMRGhostBoxFor(src);
107-
auto const lcl_box = layout.AMRToLocal(layout.AMRBoxFor(src));
112+
113+
auto const lcl_box = layout.AMRToLocal(layout.AMRBoxFor(dst));
108114

109115
if (all_primal)
110116
for (auto const lix : lcl_box)
111-
dst(lix) = src(lix);
117+
lb_view(lix) = src(lix);
112118
else
113119
for (auto const lix : lcl_box)
114-
dst(lix) = convert_to_primal(src, layout, lix, qty);
120+
lb_view(lix) = convert_to_primal(src, layout, lix, src.physicalQuantity());
115121

116122
return dst;
117123
}

0 commit comments

Comments
 (0)