Skip to content
Snippets Groups Projects
Commit e8909dcb authored by Abhinav Singh's avatar Abhinav Singh
Browse files

Fixing new Examples

parent 8c025232
No related branches found
No related tags found
1 merge request!9Fixing odeint Examples
Pipeline #3017 failed
...@@ -272,6 +272,8 @@ int main(int argc, char* argv[]) { ...@@ -272,6 +272,8 @@ int main(int argc, char* argv[]) {
// Now we initialize the grid with a filled circle. Outside the circle, the value of Phi_0 will be -1, inside +1. // Now we initialize the grid with a filled circle. Outside the circle, the value of Phi_0 will be -1, inside +1.
//Now we construct the subsets based on the subset number. //Now we construct the subsets based on the subset number.
dist_vector_subset_type Particles_bulk(Particles, 0); dist_vector_subset_type Particles_bulk(Particles, 0);
dist_vector_subset_type Particles_boundary(Particles, 1);
//We cast the global pointers to Particles and Particles_bulk as expected by the RHS functor. //We cast the global pointers to Particles and Particles_bulk as expected by the RHS functor.
PointerDistGlobal = (void *) &Particles; PointerDistGlobal = (void *) &Particles;
PointerDistSubset = (void *) &Particles_bulk; PointerDistSubset = (void *) &Particles_bulk;
...@@ -360,13 +362,15 @@ int main(int argc, char* argv[]) { ...@@ -360,13 +362,15 @@ int main(int argc, char* argv[]) {
//! @cond [OdeintT] @endcond //! @cond [OdeintT] @endcond
// Now we will create a time loop for controlling the step size ourselves but call odeint to do the 4 stages of RK4. // Now we will create a time loop for controlling the step size ourselves but call odeint to do the 4 stages of RK4.
int ctr = 0; int ctr = 0;
double t = 0, tf = 1e-1, dt = 1e-2; double t = 0, tf = 1, dt = 1e-2;
while (t < tf) { while (t < tf) {
//computing the velocity at the current position and at time step t (only in the bulk, so boundary remains 0). //computing the velocity at the current position and at time step t (only in the bulk, so boundary remains 0).
V_bulk[x] = -Pos[y] * exp(-10.0 * (Pos[x] * Pos[x] + Pos[y] * Pos[y])); V_bulk[x] = -Pos[y] * exp(-10.0 * (Pos[x] * Pos[x] + Pos[y] * Pos[y]));
V_bulk[y] = Pos[x] * exp(-10.0 * (Pos[x] * Pos[x] + Pos[y] * Pos[y])); V_bulk[y] = Pos[x] * exp(-10.0 * (Pos[x] * Pos[x] + Pos[y] * Pos[y]));
//Observing the state //Observing the state
Particles.deleteGhost();
Particles.write_frame("PDE_Sol", ctr); Particles.write_frame("PDE_Sol", ctr);
Particles.ghost_get<0>();
//calling rk4 with the function System, the state_type, current time t and the stepsize dt. It computes one step of the RK4. //calling rk4 with the function System, the state_type, current time t and the stepsize dt. It computes one step of the RK4.
Odeint_rk4.do_step(System, X, t, dt); Odeint_rk4.do_step(System, X, t, dt);
//Copying back the step, updating only the bulk values. //Copying back the step, updating only the bulk values.
...@@ -379,8 +383,12 @@ int main(int argc, char* argv[]) { ...@@ -379,8 +383,12 @@ int main(int argc, char* argv[]) {
Particles.ghost_get<0>(); Particles.ghost_get<0>();
//We update the subset and operators as the particles moved. //We update the subset and operators as the particles moved.
Particles_bulk.update(); Particles_bulk.update();
Particles_boundary.update();
Dxx.update(Particles); Dxx.update(Particles);
Dyy.update(Particles); Dyy.update(Particles);
//Reinitialzing as distributed size can change.
X.data.get<0>()=C[x];
X.data.get<1>()=C[y];
ctr++; ctr++;
t += dt; t += dt;
} //time loop end } //time loop end
...@@ -398,9 +406,4 @@ int main(int argc, char* argv[]) { ...@@ -398,9 +406,4 @@ int main(int argc, char* argv[]) {
* ## Full code ## {#odeint_c_full} * ## Full code ## {#odeint_c_full}
* *
* @include example/Numerics/OdeInt/main.cpp * @include example/Numerics/OdeInt/main.cpp
*/ */
\ No newline at end of file
...@@ -74,9 +74,9 @@ ...@@ -74,9 +74,9 @@
constexpr int x = 0; constexpr int x = 0;
constexpr int y = 1; constexpr int y = 1;
double dt; double dt=1e-2,tf=1.0;
void *PointerDistGlobal, *PointerDistSubset; void *PointerDistGlobal, *PointerDistSubset,*PointerDistSubset2;
typedef aggregate<VectorS<2, double>, VectorS<2, double>, VectorS<2, double>> Property_type; typedef aggregate<VectorS<2, double>, VectorS<2, double>, VectorS<2, double>> Property_type;
typedef vector_dist_ws<2, double, Property_type> dist_vector_type; typedef vector_dist_ws<2, double, Property_type> dist_vector_type;
...@@ -203,6 +203,7 @@ struct ObserverFunctor { ...@@ -203,6 +203,7 @@ struct ObserverFunctor {
//Casting the pointers to OpenFPM vector distributions //Casting the pointers to OpenFPM vector distributions
dist_vector_type &Particles = *(dist_vector_type *) PointerDistGlobal; dist_vector_type &Particles = *(dist_vector_type *) PointerDistGlobal;
dist_vector_subset_type &Particles_bulk = *(dist_vector_subset_type *) PointerDistSubset; dist_vector_subset_type &Particles_bulk = *(dist_vector_subset_type *) PointerDistSubset;
dist_vector_subset_type &Particles_boundary = *(dist_vector_subset_type *) PointerDistSubset2;
//Aliasing the position and properties. //Aliasing the position and properties.
auto Pos = getV<PROP_POS>(Particles); auto Pos = getV<PROP_POS>(Particles);
...@@ -225,6 +226,7 @@ struct ObserverFunctor { ...@@ -225,6 +226,7 @@ struct ObserverFunctor {
Particles.ghost_get<0>(); Particles.ghost_get<0>();
//Updating the subset and operators based on new positions //Updating the subset and operators based on new positions
Particles_bulk.update(); Particles_bulk.update();
Particles_boundary.update();
Dxx.update(Particles); Dxx.update(Particles);
Dyy.update(Particles); Dyy.update(Particles);
...@@ -232,9 +234,11 @@ struct ObserverFunctor { ...@@ -232,9 +234,11 @@ struct ObserverFunctor {
X.data.get<0>() = Concentration[x]; X.data.get<0>() = Concentration[x];
X.data.get<1>() = Concentration[y]; X.data.get<1>() = Concentration[y];
//counting the step number //counting the step number
ctr++;
}
}
ctr++;
std::cout<<"Taking a step at t="<<t<<" with dt="<<t-t_old<<std::endl;
t_old=t;
//Writing the data after every step //Writing the data after every step
Particles.deleteGhost(); Particles.deleteGhost();
Particles.write_frame("PDE_sol", ctr); Particles.write_frame("PDE_sol", ctr);
...@@ -335,9 +339,12 @@ int main(int argc, char *argv[]) ...@@ -335,9 +339,12 @@ int main(int argc, char *argv[])
// Now we initialize the grid with a filled circle. Outside the circle, the value of Phi_0 will be -1, inside +1. // Now we initialize the grid with a filled circle. Outside the circle, the value of Phi_0 will be -1, inside +1.
//Now we construct the subsets based on the subset number. //Now we construct the subsets based on the subset number.
dist_vector_subset_type Particles_bulk(Particles, 0); dist_vector_subset_type Particles_bulk(Particles, 0);
dist_vector_subset_type Particles_boundary(Particles, 1);
//We cast the global pointers to Particles and Particles_bulk as expected by the RHS functor. //We cast the global pointers to Particles and Particles_bulk as expected by the RHS functor.
PointerDistGlobal = (void *) &Particles; PointerDistGlobal = (void *) &Particles;
PointerDistSubset = (void *) &Particles_bulk; PointerDistSubset = (void *) &Particles_bulk;
PointerDistSubset2 = (void *) &Particles_boundary;
//! @cond [Pointer2Init] @endcond //! @cond [Pointer2Init] @endcond
/** /**
...@@ -419,7 +426,6 @@ int main(int argc, char *argv[]) ...@@ -419,7 +426,6 @@ int main(int argc, char *argv[])
*/ */
//////////////////////////////////////////////////////////////////////////////////////////////////////////////////// ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
//! @cond [OdeintTCall] @endcond //! @cond [OdeintTCall] @endcond
double t = 0, tf = 1e-1, dt = 1e-2;
std::vector<double> inter_times; // vector to store intermediate time steps taken by odeint. std::vector<double> inter_times; // vector to store intermediate time steps taken by odeint.
size_t steps = boost::numeric::odeint::integrate_const(Odeint_rk4, System, X, 0.0, tf, dt, ObserveAndUpdate); size_t steps = boost::numeric::odeint::integrate_const(Odeint_rk4, System, X, 0.0, tf, dt, ObserveAndUpdate);
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment