Hi all,
I have to solve the continuity equation in 1D using Implicit Runge Kuta
method.
The formula I made is this: d(U_t+1)/dt = [-grad(\Gamma_e)+R_i]/exp(U_t)
\Gamma_e ,R_i : Load the value calculated by another function.
U_t : Solution calculated at the previous time.
An error occurs here. When this function is repeated, the value I load
changes.
U_t, the value that should be divided, becomes 0 and an error occurs.
I do not know why this value changes. I just want to get a fixed value.
template<int dim>
Vector<double> Electron<dim>::evaluate_diffusion_energy (const Vector<double>
&y)
{
Vector<double> tmp(dof_handler.n_dofs());
tmp = 0.;
system_matrix.vmult(tmp,y); //tmp=(A-B-C)*y
const QGauss<dim> quadrature_formula(fe_degree + 1);
FEValues<dim> fe_values(fe,
quadrature_formula,
update_values |update_gradients|
update_quadrature_points |
update_JxW_values);
const unsigned int dofs_per_cell = fe.dofs_per_cell;
const unsigned int n_q_points = quadrature_formula.size();
std::vector<Tensor<1, dim>> grad_V(n_q_points);
std::vector<Tensor<1, dim>> grad_energyflux(n_q_points);
Vector<double> cell_source(dofs_per_cell);
std::vector<types::global_dof_index> local_dof_indices(dofs_per_cell);
for (const auto &cell : dof_handler.active_cell_iterators())
{
cell->get_dof_indices(local_dof_indices);
cell_source = 0.;
fe_values.reinit(cell);
fe_values.get_function_gradients(Potential,grad_V);
fe_values.get_function_gradients(EnergyFlux,grad_energyflux);
for (unsigned int q_point = 0; q_point < n_q_points; ++q_point)
{
double E = coeffi_E(local_dof_indices,grad_V[q_point]);
for (unsigned int i = 0; i < dofs_per_cell; ++i)
cell_source(i) += fe_values.JxW(q_point)*fe_values.shape_value(i,
q_point)* // phi_i(x)
(E-grad_energyflux[i][0])/exp(initial_EnergyDensity[local_dof_indices[i]]);
}
constraints_2.distribute_local_to_global(cell_source,
local_dof_indices,
tmp); //(A-B-C)*y+E
}
Vector<double> value(dof_handler.n_dofs());
inverse_mass_matrix.vmult(value,tmp); //value=M^-1*tmp
return value;
}
I look forward to hearing from you, thank you.
Best regards,
sora
--
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see
https://groups.google.com/d/forum/dealii?hl=en
---
You received this message because you are subscribed to the Google Groups
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email
to [email protected].
To view this discussion on the web visit
https://groups.google.com/d/msgid/dealii/b4760e9d-b85e-4877-b2a3-8c4665d33472%40googlegroups.com.
