Graph Framework
Loading...
Searching...
No Matches
absorption.hpp
Go to the documentation of this file.
1//------------------------------------------------------------------------------
6//------------------------------------------------------------------------------
7//------------------------------------------------------------------------------
90//------------------------------------------------------------------------------
91#ifndef absorption_h
92#define absorption_h
93
94#include <thread>
95
96#include "newton.hpp"
97#include "output.hpp"
98#include "dispersion.hpp"
99
101namespace absorption {
102//******************************************************************************
103// Base class
104//******************************************************************************
105//------------------------------------------------------------------------------
110//------------------------------------------------------------------------------
111 template<jit::complex_scalar T, bool SAFE_MATH=true>
112 class method {
113 public:
115 typedef T base;
117 static constexpr bool safe_math = SAFE_MATH;
118
119//------------------------------------------------------------------------------
121//------------------------------------------------------------------------------
122 virtual void compile()=0;
123
124//------------------------------------------------------------------------------
128//------------------------------------------------------------------------------
129 virtual void run(const size_t time_index)=0;
130 };
131
133 template<class A>
134 concept model = std::is_base_of<method<typename A::base, A::safe_math>, A>::value;
135
136//******************************************************************************
137// Root finder.
138//******************************************************************************
139//------------------------------------------------------------------------------
144//------------------------------------------------------------------------------
145 template<jit::complex_scalar T, bool SAFE_MATH=true>
146 class root_finder : public method<T, SAFE_MATH> {
147 private:
149 graph::shared_leaf<T, SAFE_MATH> kamp;
150
152 graph::shared_leaf<T, SAFE_MATH> w;
154 graph::shared_leaf<T, SAFE_MATH> kx;
156 graph::shared_leaf<T, SAFE_MATH> ky;
158 graph::shared_leaf<T, SAFE_MATH> kz;
160 graph::shared_leaf<T, SAFE_MATH> x;
162 graph::shared_leaf<T, SAFE_MATH> y;
164 graph::shared_leaf<T, SAFE_MATH> z;
166 graph::shared_leaf<T, SAFE_MATH> t;
167
169 graph::shared_leaf<T, SAFE_MATH> residual;
170
172 workflow::manager<T, SAFE_MATH> work;
174 const size_t index;
175
177 output::result_file file;
179 output::data_set<T> dataset;
180
182 std::thread sync;
183
184 public:
185//------------------------------------------------------------------------------
200//------------------------------------------------------------------------------
201 root_finder(graph::shared_leaf<T, SAFE_MATH> kamp,
202 graph::shared_leaf<T, SAFE_MATH> w,
203 graph::shared_leaf<T, SAFE_MATH> kx,
204 graph::shared_leaf<T, SAFE_MATH> ky,
205 graph::shared_leaf<T, SAFE_MATH> kz,
206 graph::shared_leaf<T, SAFE_MATH> x,
207 graph::shared_leaf<T, SAFE_MATH> y,
208 graph::shared_leaf<T, SAFE_MATH> z,
209 graph::shared_leaf<T, SAFE_MATH> t,
210 equilibrium::shared<T, SAFE_MATH> &eq,
211 const std::string &filename="",
212 const size_t index=0) :
213 kamp(kamp), w(w), kx(kx), ky(ky), kz(kz), x(x), y(y), z(z), t(t),
214 work(index), index(index), file(filename), dataset(file), sync([]{}) {
215 auto kvec = kx*eq->get_esup1(x, y, z)
216 + ky*eq->get_esup2(x, y, z)
217 + kz*eq->get_esup3(x, y, z);
218 auto klen = kvec->length();
219
220 auto kamp_vec = kamp*kvec/klen;
221
222 graph::input_nodes<T, SAFE_MATH> inputs = {
223 graph::variable_cast(this->kamp),
224 graph::variable_cast(this->kx),
225 graph::variable_cast(this->ky),
226 graph::variable_cast(this->kz),
227 graph::variable_cast(this->x),
228 graph::variable_cast(this->y),
229 graph::variable_cast(this->z)
230 };
231
232 graph::map_nodes<T, SAFE_MATH> setters = {
233 {graph::zero<T, SAFE_MATH> (), graph::variable_cast(this->kamp)}
234 };
235
236 work.add_item(inputs, {}, setters, NULL,
237 "root_find_init_kernel", inputs.back()->size());
238
239 inputs.push_back(graph::variable_cast(this->t));
240 inputs.push_back(graph::variable_cast(this->w));
241
242 auto D = dispersion::hot_plasma<T,
243 dispersion::z_erfi<T, SAFE_MATH>,
244 SAFE_MATH>().D(w,
245 kvec + kamp_vec,
246 x, y, z, t, eq);
247
248 solver::newton(work, {kamp}, inputs, {D},
249 graph::shared_random_state<T, SAFE_MATH> ());
250
251 inputs = {
252 graph::variable_cast(this->kamp),
253 graph::variable_cast(this->kx),
254 graph::variable_cast(this->ky),
255 graph::variable_cast(this->kz),
256 graph::variable_cast(this->x),
257 graph::variable_cast(this->y),
258 graph::variable_cast(this->z)
259 };
260 setters = {
261 {klen + kamp, graph::variable_cast(this->kamp)}
262 };
263 work.add_item(inputs, {}, setters, NULL,
264 "final_kamp", inputs.back()->size());
265 }
266
267//------------------------------------------------------------------------------
269//------------------------------------------------------------------------------
270 ~root_finder() {
271 sync.join();
272 }
273
274//------------------------------------------------------------------------------
276//------------------------------------------------------------------------------
277 void compile() {
278 work.compile();
279
280 dataset.create_variable(file, "kamp", this->kamp, work.get_context());
281
282 dataset.reference_variable(file, "w", graph::variable_cast(this->w));
283 dataset.reference_variable(file, "kx", graph::variable_cast(this->kx));
284 dataset.reference_variable(file, "ky", graph::variable_cast(this->ky));
285 dataset.reference_variable(file, "kz", graph::variable_cast(this->kz));
286 dataset.reference_variable(file, "x", graph::variable_cast(this->x));
287 dataset.reference_variable(file, "y", graph::variable_cast(this->y));
288 dataset.reference_variable(file, "z", graph::variable_cast(this->z));
289 dataset.reference_variable(file, "time", graph::variable_cast(this->t));
290 file.end_define_mode();
291 }
292
293//------------------------------------------------------------------------------
297//------------------------------------------------------------------------------
298 void run(const size_t time_index) {
299 dataset.read(file, time_index);
300 work.copy_to_device(w, graph::variable_cast(this->w)->data());
301 work.copy_to_device(kx, graph::variable_cast(this->kx)->data());
302 work.copy_to_device(ky, graph::variable_cast(this->ky)->data());
303 work.copy_to_device(kz, graph::variable_cast(this->kz)->data());
304 work.copy_to_device(x, graph::variable_cast(this->x)->data());
305 work.copy_to_device(y, graph::variable_cast(this->y)->data());
306 work.copy_to_device(z, graph::variable_cast(this->z)->data());
307 work.copy_to_device(t, graph::variable_cast(this->t)->data());
308
309 work.run();
310
311 sync.join();
312 work.wait();
313 sync = std::thread([this] (const size_t i) -> void {
314 dataset.write(file, i);
315 }, time_index);
316 }
317 };
318
319//******************************************************************************
320// Weak Damping.
321//******************************************************************************
322//------------------------------------------------------------------------------
327//------------------------------------------------------------------------------
328 template<jit::complex_scalar T, bool SAFE_MATH=true>
329 class weak_damping : public method<T, SAFE_MATH> {
330 private:
332 graph::shared_leaf<T, SAFE_MATH> kamp;
333
335 graph::shared_leaf<T, SAFE_MATH> w;
337 graph::shared_leaf<T, SAFE_MATH> kx;
339 graph::shared_leaf<T, SAFE_MATH> ky;
341 graph::shared_leaf<T, SAFE_MATH> kz;
343 graph::shared_leaf<T, SAFE_MATH> x;
345 graph::shared_leaf<T, SAFE_MATH> y;
347 graph::shared_leaf<T, SAFE_MATH> z;
349 graph::shared_leaf<T, SAFE_MATH> t;
350
352 workflow::manager<T, SAFE_MATH> work;
353
355 const size_t index;
356
358 output::result_file file;
360 output::data_set<T> dataset;
361
363 std::thread sync;
364
365 public:
366//------------------------------------------------------------------------------
381//------------------------------------------------------------------------------
382 weak_damping(graph::shared_leaf<T, SAFE_MATH> kamp,
383 graph::shared_leaf<T, SAFE_MATH> w,
384 graph::shared_leaf<T, SAFE_MATH> kx,
385 graph::shared_leaf<T, SAFE_MATH> ky,
386 graph::shared_leaf<T, SAFE_MATH> kz,
387 graph::shared_leaf<T, SAFE_MATH> x,
388 graph::shared_leaf<T, SAFE_MATH> y,
389 graph::shared_leaf<T, SAFE_MATH> z,
390 graph::shared_leaf<T, SAFE_MATH> t,
391 equilibrium::shared<T, SAFE_MATH> &eq,
392 const std::string &filename="",
393 const size_t index=0) :
394 kamp(kamp), w(w), kx(kx), ky(ky), kz(kz), x(x), y(y), z(z), t(t),
395 work(index), index(index), file(filename), dataset(file), sync([]{}) {
396 auto k_vec = kx*eq->get_esup1(x, y, z)
397 + ky*eq->get_esup2(x, y, z)
398 + kz*eq->get_esup3(x, y, z);
399 auto k_unit = k_vec->unit();
400
401 auto Dc = dispersion::cold_plasma_expansion<T, SAFE_MATH> ().D(w, k_vec,
402 x, y, z, t, eq);
403 auto Dw = dispersion::hot_plasma_expansion<T,
404 dispersion::z_erfi<T, SAFE_MATH>,
405 SAFE_MATH> ().D(w, k_vec,
406 x, y, z, t, eq);
407
408 auto kamp1 = k_vec->length()
409 - Dw/k_unit->dot(Dc->df(kx)*eq->get_esup1(x, y, z) +
410 Dc->df(ky)*eq->get_esup2(x, y, z) +
411 Dc->df(kz)*eq->get_esup3(x, y, z));
412
413 graph::input_nodes<T, SAFE_MATH> inputs = {
414 graph::variable_cast(this->kamp),
415 graph::variable_cast(this->kx),
416 graph::variable_cast(this->ky),
417 graph::variable_cast(this->kz),
418 graph::variable_cast(this->x),
419 graph::variable_cast(this->y),
420 graph::variable_cast(this->z),
421 graph::variable_cast(this->t),
422 graph::variable_cast(this->w)
423 };
424
425 graph::map_nodes<T, SAFE_MATH> setters = {
426 {kamp1, graph::variable_cast(this->kamp)}
427 };
428
429 work.add_item(inputs, {}, setters,
430 graph::shared_random_state<T, SAFE_MATH> (),
431 "weak_damping_kimg_kernel", inputs.back()->size());
432 }
433
434//------------------------------------------------------------------------------
436//------------------------------------------------------------------------------
437 ~weak_damping() {
438 sync.join();
439 }
440
441//------------------------------------------------------------------------------
443//------------------------------------------------------------------------------
444 void compile() {
445 work.compile();
446
447 dataset.create_variable(file, "kamp", this->kamp, work.get_context());
448
449 dataset.reference_variable(file, "w", graph::variable_cast(this->w));
450 dataset.reference_variable(file, "kx", graph::variable_cast(this->kx));
451 dataset.reference_variable(file, "ky", graph::variable_cast(this->ky));
452 dataset.reference_variable(file, "kz", graph::variable_cast(this->kz));
453 dataset.reference_variable(file, "x", graph::variable_cast(this->x));
454 dataset.reference_variable(file, "y", graph::variable_cast(this->y));
455 dataset.reference_variable(file, "z", graph::variable_cast(this->z));
456 dataset.reference_variable(file, "time", graph::variable_cast(this->t));
457 file.end_define_mode();
458 }
459
460//------------------------------------------------------------------------------
464//------------------------------------------------------------------------------
465 void run(const size_t time_index) {
466 dataset.read(file, time_index);
467 work.copy_to_device(w, graph::variable_cast(this->w)->data());
468 work.copy_to_device(kx, graph::variable_cast(this->kx)->data());
469 work.copy_to_device(ky, graph::variable_cast(this->ky)->data());
470 work.copy_to_device(kz, graph::variable_cast(this->kz)->data());
471 work.copy_to_device(x, graph::variable_cast(this->x)->data());
472 work.copy_to_device(y, graph::variable_cast(this->y)->data());
473 work.copy_to_device(z, graph::variable_cast(this->z)->data());
474 work.copy_to_device(t, graph::variable_cast(this->t)->data());
475
476 work.run();
477
478 sync.join();
479 work.wait();
480 sync = std::thread([this] (const size_t i) -> void {
481 dataset.write(file, i);
482 }, time_index);
483 }
484 };
485}
486
487#endif /* absorption_h */
absorption::method
Base class for absorption models.
Definition absorption.hpp:112
absorption::method::safe_math
static constexpr bool safe_math
Retrieve template parameter of safe math.
Definition absorption.hpp:117
absorption::method::base
T base
Type def to retrieve the backend base type.
Definition absorption.hpp:115
absorption::method::run
virtual void run(const size_t time_index)=0
Run the workflow.
absorption::method::compile
virtual void compile()=0
Compile the work items.
absorption::root_finder
Class interface for the root finder.
Definition absorption.hpp:146
absorption::root_finder::root_finder
root_finder(graph::shared_leaf< T, SAFE_MATH > kamp, graph::shared_leaf< T, SAFE_MATH > w, graph::shared_leaf< T, SAFE_MATH > kx, graph::shared_leaf< T, SAFE_MATH > ky, graph::shared_leaf< T, SAFE_MATH > kz, graph::shared_leaf< T, SAFE_MATH > x, graph::shared_leaf< T, SAFE_MATH > y, graph::shared_leaf< T, SAFE_MATH > z, graph::shared_leaf< T, SAFE_MATH > t, equilibrium::shared< T, SAFE_MATH > &eq, const std::string &filename="", const size_t index=0)
Constructor for root finding.
Definition absorption.hpp:201
absorption::root_finder::run
void run(const size_t time_index)
Run the workflow.
Definition absorption.hpp:298
absorption::root_finder::compile
void compile()
Compile the work items.
Definition absorption.hpp:277
absorption::root_finder::~root_finder
~root_finder()
Destructor.
Definition absorption.hpp:270
absorption::weak_damping
Class interface weak damping approximation.
Definition absorption.hpp:329
absorption::weak_damping::run
void run(const size_t time_index)
Run the workflow.
Definition absorption.hpp:465
absorption::weak_damping::~weak_damping
~weak_damping()
Destructor.
Definition absorption.hpp:437
absorption::weak_damping::compile
void compile()
Compile the work items.
Definition absorption.hpp:444
absorption::weak_damping::weak_damping
weak_damping(graph::shared_leaf< T, SAFE_MATH > kamp, graph::shared_leaf< T, SAFE_MATH > w, graph::shared_leaf< T, SAFE_MATH > kx, graph::shared_leaf< T, SAFE_MATH > ky, graph::shared_leaf< T, SAFE_MATH > kz, graph::shared_leaf< T, SAFE_MATH > x, graph::shared_leaf< T, SAFE_MATH > y, graph::shared_leaf< T, SAFE_MATH > z, graph::shared_leaf< T, SAFE_MATH > t, equilibrium::shared< T, SAFE_MATH > &eq, const std::string &filename="", const size_t index=0)
Constructor for weak damping.
Definition absorption.hpp:382
dispersion::cold_plasma_expansion
Cold Plasma expansion dispersion function.
Definition dispersion.hpp:1018
dispersion::cold_plasma_expansion::D
virtual graph::shared_leaf< T, SAFE_MATH > D(graph::shared_leaf< T, SAFE_MATH > w, graph::shared_vector< T, SAFE_MATH > k_vec, graph::shared_leaf< T, SAFE_MATH > x, graph::shared_leaf< T, SAFE_MATH > y, graph::shared_leaf< T, SAFE_MATH > z, graph::shared_leaf< T, SAFE_MATH > t, equilibrium::shared< T, SAFE_MATH > &eq)
Cold Plasma expansion Dispersion function.
Definition dispersion.hpp:1043
dispersion::hot_plasma_expansion
Hot Plasma Expansion Dispersion function.
Definition dispersion.hpp:1209
dispersion::hot_plasma
Hot Plasma Dispersion function.
Definition dispersion.hpp:1100
dispersion::z_erfi
Class interface to build dispersion relation functions.
Definition dispersion.hpp:289
output::data_set
Class representing a netcdf dataset.
Definition output.hpp:166
output::data_set::read
void read(const result_file &result, const size_t index)
Read step.
Definition output.hpp:412
output::data_set::reference_variable
void reference_variable(const result_file &result, const std::string &name, graph::shared_variable< T, SAFE_MATH > &&node)
Load reference.
Definition output.hpp:285
output::data_set::create_variable
void create_variable(const result_file &result, const std::string &name, graph::shared_leaf< T, SAFE_MATH > &node, jit::context< T, SAFE_MATH > &context)
Create a variable.
Definition output.hpp:260
output::data_set::write
void write(const result_file &result)
Write step.
Definition output.hpp:353
output::result_file
Class representing a netcdf based output file.
Definition output.hpp:32
output::result_file::end_define_mode
void end_define_mode() const
End define mode.
Definition output.hpp:94
workflow::manager
Class representing a workflow manager.
Definition workflow.hpp:215
workflow::manager::run
void run()
Run work items.
Definition workflow.hpp:359
workflow::manager::wait
void wait()
Wait for GPU queue to finish.
Definition workflow.hpp:368
workflow::manager::add_item
void add_item(graph::input_nodes< T, SAFE_MATH > in, graph::output_nodes< T, SAFE_MATH > out, graph::map_nodes< T, SAFE_MATH > maps, graph::shared_random_state< T, SAFE_MATH > state, const std::string name, const size_t size)
Add a workflow item.
Definition workflow.hpp:271
workflow::manager::copy_to_device
void copy_to_device(graph::shared_leaf< T, SAFE_MATH > &node, T *destination)
Copy buffer contents to the device.
Definition workflow.hpp:378
workflow::manager::get_context
jit::context< T, SAFE_MATH > & get_context()
Get the jit context.
Definition workflow.hpp:422
workflow::manager::compile
void compile()
Compile the workflow items.
Definition workflow.hpp:336
absorption::model
Solver method concept.
Definition absorption.hpp:134
dispersion.hpp
Base class for a dispersion relation.
absorption
Namespace for power absorption models.
Definition absorption.hpp:101
equilibrium::shared
std::shared_ptr< generic< T, SAFE_MATH > > shared
Convenience type alias for shared equilibria.
Definition equilibrium.hpp:470
graph::zero
constexpr shared_leaf< T, SAFE_MATH > zero()
Forward declare for zero.
Definition node.hpp:986
graph::input_nodes
std::vector< shared_variable< T, SAFE_MATH > > input_nodes
Convenience type alias for a vector of inputs.
Definition node.hpp:1711
graph::variable_cast
shared_variable< T, SAFE_MATH > variable_cast(shared_leaf< T, SAFE_MATH > x)
Cast to a variable node.
Definition node.hpp:1727
graph::shared_leaf
std::shared_ptr< leaf_node< T, SAFE_MATH > > shared_leaf
Convenience type alias for shared leaf nodes.
Definition node.hpp:676
graph::map_nodes
std::vector< std::pair< shared_leaf< T, SAFE_MATH >, shared_variable< T, SAFE_MATH > > > map_nodes
Convenience type alias for mapping end codes back to inputs.
Definition node.hpp:1715
solver::newton
void newton(workflow::manager< T, SAFE_MATH > &work, graph::output_nodes< T, SAFE_MATH > vars, graph::input_nodes< T, SAFE_MATH > inputs, graph::shared_leaf< T, SAFE_MATH > func, graph::shared_random_state< T, SAFE_MATH > state, const T tolerance=1.0E-30, const size_t max_iterations=1000, const T step=1.0)
Determine the value of vars to minimize the loss function.
Definition newton.hpp:34
newton.hpp
Sets up the kernel for a newtons method.
output.hpp
Implements output files in a netcdf format.