[deal.II] Re: increment Dirichlet boundary condition

[Bryukhanov Ilya]

I tried to use code from step-22 and from the sources you recommended.
The problem is that the solution is zero over the body except the boundary 
where
Dirichlet conditions were set. It is just one iteration. I don't know why 
they are not in use.
Briefly, my code consists of the following steps which are perfomed in turn 
(from 1 to 5)
1) Setup variables and constraints
dof_handler.distribute_dofs (fe);
{
constraints.clear ();
DoFTools::make_hanging_node_constraints(dof_handler, constraints);
VectorTools::interpolate_boundary_values(dof_handler,
 1,
 boundary_u_handler[1],
 constraints);

VectorTools::interpolate_boundary_values(dof_handler,
 3,
 boundary_u_handler[3],
 constraints);
}
constraints.close ();
DynamicSparsityPattern dsp(dof_handler.n_dofs(), dof_handler.n_dofs());
DoFTools::make_sparsity_pattern(dof_handler, dsp, constraints, false);
sparsity_pattern.copy_from (dsp);
system_matrix.reinit (sparsity_pattern);
solution.reinit (dof_handler.n_dofs());
system_rhs.reinit (dof_handler.n_dofs()); 

2) Assembling FEM matrix
FEValues 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();

FullMatrix   cell_matrix (dofs_per_cell, dofs_per_cell);
Vector   cell_rhs (dofs_per_cell);
std::vector local_dof_indices (dofs_per_cell);

// Loop over all cells:
for (auto cell = dof_handler.begin_active(); cell!=dof_handler.end(); 
++cell) {
cell_matrix = 0;
cell_rhs = 0;
fe_values.reinit(cell);
for (unsigned int i = 0; i < dofs_per_cell; ++i) {
for (unsigned int j = 0; j < dofs_per_cell; ++j) {
for (unsigned int q_point = 0; q_point < n_q_points;
 ++q_point) {
const SymmetricTensor<2, dim>
eps_phi_i = get_strain(fe_values, i, q_point),
eps_phi_j = get_strain(fe_values, j, q_point);
cell_matrix(i, j) +=
(eps_phi_i * stress_strain_tensor * eps_phi_j
 *
 fe_values.JxW(q_point));
}
}
}
cell->get_dof_indices (local_dof_indices);
constraints.
distribute_local_to_global(cell_matrix, cell_rhs,
   local_dof_indices,
   system_matrix, system_rhs);
}

3) Assembling RHS vector
 FEValues fe_values (fe, quadrature_formula,
 update_values   | update_gradients |
 update_quadrature_points | update_JxW_values);

FEFaceValues fe_face_values (fe, face_quadrature_formula,
  update_values   | 
update_normal_vectors |
  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();
const unsigned int   n_face_q_points = face_quadrature_formula.size();
FullMatrix   cell_matrix (dofs_per_cell, dofs_per_cell);

// right hand side vector for the cell
Vector   cell_rhs (dofs_per_cell);
std::vector local_dof_indices (dofs_per_cell);
std::vector rhs_values (n_q_points,
 Vector(dim));
// reinit right hand side
system_rhs.reinit (dof_handler.n_dofs());
// loop over all cells:
for (auto cell = dof_handler.begin_active(); cell!=dof_handler.end(); 
++cell) {
cell_matrix = 0;
cell_rhs = 0;
fe_values.reinit (cell);
// Assembling the right hand side boundary terms
for (unsigned int face_number=0;
 face_number neumann_vector_value;
if (cell->face(face_number)->at_boundary()) {
int id = cell->face(face_number)->boundary_id();
if (boundary_force.find(id) != boundary_force.end()) {
fe_face_values.reinit(cell, face_number);
// quadrature
for (unsigned int i = 0; i < dofs_per_cell; ++i) {
const unsigned int component_i = 
fe.system_to_component_index(i).first;
for (unsigned int q = 0; q < n_face_q_points; ++q) {
  

[deal.II] Re: increment Dirichlet boundary condition

[Bryukhanov Ilya]

I tried to use code from step-22 and from the sources you recommended.
The problem is that the solution is zero over the body except the boundary 
where
Dirichlet conditions were set. It is just one iteration. I don't know why 
they are not in use.
Briefly, my code consists of the following steps which are perfomed in turn 
(from 1 to 5)
1) Setup variables and constraints
dof_handler.distribute_dofs (fe);
{
constraints.clear ();
DoFTools::make_hanging_node_constraints(dof_handler, constraints);
VectorTools::interpolate_boundary_values(dof_handler,
 1,
 boundary_u_handler[1],
 constraints);

VectorTools::interpolate_boundary_values(dof_handler,
 3,
 boundary_u_handler[3],
 constraints);
}
constraints.close ();
DynamicSparsityPattern dsp(dof_handler.n_dofs(), dof_handler.n_dofs());
DoFTools::make_sparsity_pattern(dof_handler, dsp, constraints, false);
sparsity_pattern.copy_from (dsp);
system_matrix.reinit (sparsity_pattern);
solution.reinit (dof_handler.n_dofs());
system_rhs.reinit (dof_handler.n_dofs()); 

2) Assembling FEM matrix
FEValues 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();

FullMatrix   cell_matrix (dofs_per_cell, dofs_per_cell);
Vector   cell_rhs (dofs_per_cell);
std::vector local_dof_indices (dofs_per_cell);

// Loop over all cells:
for (auto cell = dof_handler.begin_active(); cell!=dof_handler.end(); 
++cell) {
cell_matrix = 0;
cell_rhs = 0;
fe_values.reinit(cell);
for (unsigned int i = 0; i < dofs_per_cell; ++i) {
for (unsigned int j = 0; j < dofs_per_cell; ++j) {
for (unsigned int q_point = 0; q_point < n_q_points;
 ++q_point) {
const SymmetricTensor<2, dim>
eps_phi_i = get_strain(fe_values, i, q_point),
eps_phi_j = get_strain(fe_values, j, q_point);
cell_matrix(i, j) +=
(eps_phi_i * stress_strain_tensor * eps_phi_j
 *
 fe_values.JxW(q_point));
}
}
}
cell->get_dof_indices (local_dof_indices);
constraints.
distribute_local_to_global(cell_matrix, cell_rhs,
   local_dof_indices,
   system_matrix, system_rhs);
}

3) Assembling RHS vector
 FEValues fe_values (fe, quadrature_formula,
 update_values   | update_gradients |
 update_quadrature_points | update_JxW_values);

FEFaceValues fe_face_values (fe, face_quadrature_formula,
  update_values   | 
update_normal_vectors |
  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();
const unsigned int   n_face_q_points = face_quadrature_formula.size();
FullMatrix   cell_matrix (dofs_per_cell, dofs_per_cell);

// right hand side vector for the cell
Vector   cell_rhs (dofs_per_cell);
std::vector local_dof_indices (dofs_per_cell);
std::vector rhs_values (n_q_points,
 Vector(dim));
// reinit right hand side
system_rhs.reinit (dof_handler.n_dofs());
// loop over all cells:
for (auto cell = dof_handler.begin_active(); cell!=dof_handler.end(); 
++cell) {
cell_matrix = 0;
cell_rhs = 0;
fe_values.reinit (cell);
// Assembling the right hand side boundary terms
for (unsigned int face_number=0;
 face_number neumann_vector_value;
if (cell->face(face_number)->at_boundary()) {
int id = cell->face(face_number)->boundary_id();
if (boundary_force.find(id) != boundary_force.end()) {
fe_face_values.reinit(cell, face_number);
// quadrature
for (unsigned int i = 0; i < dofs_per_cell; ++i) {
const unsigned int component_i = 
fe.system_to_component_index(i).first;
for (unsigned int q = 0; q < n_face_q_points; ++q) {
  

[deal.II] Re: increment Dirichlet boundary condition

[Denis Davydov]

Hi IIiya,

On Thursday, February 22, 2018 at 11:55:33 AM UTC+1, Bryukhanov Ilya wrote:
>
> Denis, thanks a lot for the answer!
>
> I think that on the first step I can use usual 
> "interpolate_boundary_values + apply_boundary_values" functions,
>

the outline I gave above is agnostic to "first step"/"other steps", you 
would call VectorTools::interpolate_boundary_values() regardless 
stage to get constraints into the ConstraintMatrix object. If you do 
Newton-Raphson and solve for increments, at the first non-linear step you 
should have all constraints (hanging nodes + non-zero Dirichlet for a given 
step in quasi-static(?) loading), and for consequent steps only hanging 
nodes plus zero Dirichlet. 

Please read this module:
https://www.dealii.org/developer/doxygen/deal.II/group__constraints.html

This code-gallery may help in figuring out how to do 
this: 
https://www.dealii.org/developer/doxygen/deal.II/code_gallery_Quasi_static_Finite_strain_Compressible_Elasticity.html
 

 

> that transform the global matrix: columns of the constrained nodes become 
> zero columns with diagonal non-zero entries.
>
> In a time loop I have to account boundary condition without modifing 
> global matrix. While looping over dirichlet constrained cells
> I call the function  
> "hanging_node_constraints.distribute_local_to_global(cell_rhs, 
> local_dof_indices, system_rhs, cell_matrix)"
> which modifies the system_rhs vector by adding columns from constrained 
> nodes, which are actually zero in global matrix.
>
 
Nope, if your "hanging_node_constraints" indeed contains only constraints 
from hanging nodes, you won't see any contribution to the RHS
from it as constraints are homogeneous (ie. you wont' have  u_{123} = 456 
from hanging nodes).
 

>  
> Also I need to add components to rhs vector from the force boundary 
> conditions.
>

that's a different story, just do integration over faces and distribute as 
usual.
 

>
> However, I have to somehow account Dirichlet boundary conditions as in 
> previous step I didn't use any boundary values. The function
>

use ConstraintMatrix which has non-zero Dirichlet BC constraints.
 

> "interpolate_boundary_values(dof_handler, boundary_id, DirichletHandler, 
> boundary_values)"
> gives me "boundary_values" variables that maps node numbers to dirichle 
> boundary. 
>
> But I don't know how should I use it. What should I do with that variable 
> to account boundary conditions together with "distribute_local_to_global"?
>
> Thanks a lot in advance. Sorry for my misunderstanding.
>

Regards,
Denis.

-- 
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 dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: increment Dirichlet boundary condition

[Bryukhanov Ilya]

Denis, thanks a lot for the answer!

I think that on the first step I can use usual "interpolate_boundary_values 
+ apply_boundary_values" functions,
that transform the global matrix: columns of the constrained nodes become 
zero columns with diagonal non-zero entries.

In a time loop I have to account boundary condition without modifing global 
matrix. While looping over dirichlet constrained cells
I call the function  
"hanging_node_constraints.distribute_local_to_global(cell_rhs, 
local_dof_indices, system_rhs, cell_matrix)"
which modifies the system_rhs vector by adding columns from constrained 
nodes, which are actually zero in global matrix. 
Also I need to add components to rhs vector from the force boundary 
conditions.

However, I have to somehow account Dirichlet boundary conditions as in 
previous step I didn't use any boundary values. The function
"interpolate_boundary_values(dof_handler, boundary_id, DirichletHandler, 
boundary_values)"
gives me "boundary_values" variables that maps node numbers to dirichle 
boundary. 

But I don't know how should I use it. What should I do with that variable 
to account boundary conditions together with "distribute_local_to_global"?

Thanks a lot in advance. Sorry for my misunderstanding.

On Thursday, February 22, 2018 at 1:14:41 AM UTC+3, Denis Davydov wrote:
>
> Hi Iliya,
>
> there are numerous things that are wrong
>
> On Wednesday, February 21, 2018 at 5:32:55 PM UTC+1, Bryukhanov Ilya wrote:
>>
>> Denis, thanks a lot for your answer and links that you provided. 
>>
>> However, my solution becomes in a way that the displacement values on 
>> that boundary abnormally grow
>> even without assigning new boundary condition in a time loop. I attach 
>> figures on the zero step (correct) and on the first step (wrong).
>>
>> I do the following these things on each iteration: 
>> 1) recompute rhs vector - system_rhs (which is zero for my problem),
>> 2) compute matrix terms that are on the dirichlet boundary:
>> for (auto cell = dof_handler.begin_active(); cell!=dof_handler.end(); 
>> ++cell) {
>> cell_matrix = 0; cell_rhs = 0; fe_values.reinit(cell);
>> for (unsigned int face_number=0; 
>> face_number> if (cell->face(face_number)->at_boundary()) {
>> // boundary index
>> int id = cell->face(face_number)->boundary_id();
>> if (boundary_displacement.find(id) != 
>> boundary_displacement.end()) {
>> //compute cell_matrix
>>
>> Then call function
>> hanging_node_constraints.distribute_local_to_global(cell_matrix, 
>> cell_rhs, local_dof_indices, system_matrix, system_rhs);
>>
>
> ^^ that should be hanging_node_constraints + Dirichlet 
>
>>
>> 3) Then  apply functions to account Dirichlet boundary condition:  
>>  VectorTools::interpolate_boundary_values(dof_handler, boundary_id, 
>> DirichletHandler, boundary_values);
>>
> ^^ that you would use during your setup
>  
>
>> MatrixTools::apply_boundary_values (boundary_values,
>> system_matrix,
>> solution,
>> system_rhs);
>>
>
> ^^ this is already what constraints.distribute_local_to_global() does, so 
> no need to do that.
>
> Regards
> Denis.
>  
>
>>
>>
>> Can you give some advice how to make it correct?
>>
>>
>> 
>>  
>> 
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>> On Tuesday, February 20, 2018 at 10:51:17 PM UTC+3, Denis Davydov wrote:
>>>
>>> Hi Ilya,
>>>
>>> I am not away of principally simple solution, the RHS simply contains 
>>> your bilinear form times interpolated values at constrained DoFs. 
>>> What you can do is to loop over a subset of cells (those which are at 
>>> the boundary), still assemble the local matrix and do 
>>> constraints.distribute_local_to_global() which takes global matrix and 
>>> RHS.
>>>
>>> Or you can use constrained_linear_operator 
>>> ,
>>>  
>>> see this 
>>>  forum 
>>> discussion.
>>>
>>> Please also study this 
>>>  
>>> thread, which covers similar question.
>>>
>>> Regards,
>>> Denis.
>>>
>>> On Tuesday, February 20, 2018 at 4:46:19 PM UTC+1, Bryukhanov Ilya wrote:

 Hi,

 I consider elastic problem and forces and displacements on boundaries 
 are changing with time. I don't want to compute FEM matrix
 on each step. How to correctly simulate the task? I know that there is 
 a 

[deal.II] Re: increment Dirichlet boundary condition

[Bryukhanov Ilya]

Denis, thanks a lot for the answer!

I think that on the first step I can use usual "interpolate_boundary_values 
+ apply_boundary_values" functions,
that transform the global matrix: columns of the constrained nodes become 
zero columns with diagonal non-zero entries.

In a time loop I have to account boundary condition without modifing global 
matrix. While looping over dirichlet constrained cells
I call the function  
"hanging_node_constraints.distribute_local_to_global(cell_rhs, 
local_dof_indices, system_rhs, cell_matrix)"
which modifies the system_rhs vector by adding columns from constrained 
nodes, which are actually zero in global matrix. 

However, I have to somehow account boundary conditions. The function
"interpolate_boundary_values(dof_handler, boundary_id, DirichletHandler, 
boundary_values)"
gives me "boundary_values" variables that maps node numbers to dirichle 
boundary. 
But I don't know how should I use it. What should I do with that variable 
to account boundary conditions together with "distribute_local_to_global"?

Thanks a lot in advance. Sorry for my misundersting.

On Thursday, February 22, 2018 at 1:14:41 AM UTC+3, Denis Davydov wrote:
>
> Hi Iliya,
>
> there are numerous things that are wrong
>
> On Wednesday, February 21, 2018 at 5:32:55 PM UTC+1, Bryukhanov Ilya wrote:
>>
>> Denis, thanks a lot for your answer and links that you provided. 
>>
>> However, my solution becomes in a way that the displacement values on 
>> that boundary abnormally grow
>> even without assigning new boundary condition in a time loop. I attach 
>> figures on the zero step (correct) and on the first step (wrong).
>>
>> I do the following these things on each iteration: 
>> 1) recompute rhs vector - system_rhs (which is zero for my problem),
>> 2) compute matrix terms that are on the dirichlet boundary:
>> for (auto cell = dof_handler.begin_active(); cell!=dof_handler.end(); 
>> ++cell) {
>> cell_matrix = 0; cell_rhs = 0; fe_values.reinit(cell);
>> for (unsigned int face_number=0; 
>> face_number> if (cell->face(face_number)->at_boundary()) {
>> // boundary index
>> int id = cell->face(face_number)->boundary_id();
>> if (boundary_displacement.find(id) != 
>> boundary_displacement.end()) {
>> //compute cell_matrix
>>
>> Then call function
>> hanging_node_constraints.distribute_local_to_global(cell_matrix, 
>> cell_rhs, local_dof_indices, system_matrix, system_rhs);
>>
>
> ^^ that should be hanging_node_constraints + Dirichlet 
>
>>
>> 3) Then  apply functions to account Dirichlet boundary condition:  
>>  VectorTools::interpolate_boundary_values(dof_handler, boundary_id, 
>> DirichletHandler, boundary_values);
>>
> ^^ that you would use during your setup
>  
>
>> MatrixTools::apply_boundary_values (boundary_values,
>> system_matrix,
>> solution,
>> system_rhs);
>>
>
> ^^ this is already what constraints.distribute_local_to_global() does, so 
> no need to do that.
>
> Regards
> Denis.
>  
>
>>
>>
>> Can you give some advice how to make it correct?
>>
>>
>> 
>>  
>> 
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>> On Tuesday, February 20, 2018 at 10:51:17 PM UTC+3, Denis Davydov wrote:
>>>
>>> Hi Ilya,
>>>
>>> I am not away of principally simple solution, the RHS simply contains 
>>> your bilinear form times interpolated values at constrained DoFs. 
>>> What you can do is to loop over a subset of cells (those which are at 
>>> the boundary), still assemble the local matrix and do 
>>> constraints.distribute_local_to_global() which takes global matrix and 
>>> RHS.
>>>
>>> Or you can use constrained_linear_operator 
>>> ,
>>>  
>>> see this 
>>>  forum 
>>> discussion.
>>>
>>> Please also study this 
>>>  
>>> thread, which covers similar question.
>>>
>>> Regards,
>>> Denis.
>>>
>>> On Tuesday, February 20, 2018 at 4:46:19 PM UTC+1, Bryukhanov Ilya wrote:

 Hi,

 I consider elastic problem and forces and displacements on boundaries 
 are changing with time. I don't want to compute FEM matrix
 on each step. How to correctly simulate the task? I know that there is 
 a simple solution, but my solution is incorrect :)
 I have two variables for the system_matrix and rhs_vector. I try to 
 recompute rhs_vector on 

[deal.II] Re: increment Dirichlet boundary condition

[Bryukhanov Ilya]

Denis, thanks a lot for the answer! I think that on the first step I can 
use usual "interpolate_boundary_values + apply_boundary_values" functions,
and this procedures transform the global matrix so that the columns of 
constrained nodes become zero columns with diagonal non-zero entries.
In a time loop I have to account boundary condition so the calling of 
"hanging_node_constraints.distribute_local_to_global(cell_rhs, 
local_dof_indices, system_rhs, cell_matrix)"
while looping over dirichlet constrained cells modify the system_rhs vector 
by adding columns from constrained nodes,
which are actually zero in global matrix. However, I have to somehow 
account boundary conditions. The function
"interpolate_boundary_values(dof_handler, boundary_id, DirichletHandler, 
boundary_values)"
gives me "boundary_values" variables that maps node numbers to dirichle 
boundary. But I don't know how should I use it. What should I 
do with the variable "boundary_values" to account boundary conditions in 
"distribute_local_to_global"?

Thanks a lot in advance. Sorry for my misundersting.

On Thursday, February 22, 2018 at 1:14:41 AM UTC+3, Denis Davydov wrote:
>
> Hi Iliya,
>
> there are numerous things that are wrong
>
> On Wednesday, February 21, 2018 at 5:32:55 PM UTC+1, Bryukhanov Ilya wrote:
>>
>> Denis, thanks a lot for your answer and links that you provided. 
>>
>> However, my solution becomes in a way that the displacement values on 
>> that boundary abnormally grow
>> even without assigning new boundary condition in a time loop. I attach 
>> figures on the zero step (correct) and on the first step (wrong).
>>
>> I do the following these things on each iteration: 
>> 1) recompute rhs vector - system_rhs (which is zero for my problem),
>> 2) compute matrix terms that are on the dirichlet boundary:
>> for (auto cell = dof_handler.begin_active(); cell!=dof_handler.end(); 
>> ++cell) {
>> cell_matrix = 0; cell_rhs = 0; fe_values.reinit(cell);
>> for (unsigned int face_number=0; 
>> face_number> if (cell->face(face_number)->at_boundary()) {
>> // boundary index
>> int id = cell->face(face_number)->boundary_id();
>> if (boundary_displacement.find(id) != 
>> boundary_displacement.end()) {
>> //compute cell_matrix
>>
>> Then call function
>> hanging_node_constraints.distribute_local_to_global(cell_matrix, 
>> cell_rhs, local_dof_indices, system_matrix, system_rhs);
>>
>
> ^^ that should be hanging_node_constraints + Dirichlet 
>
>>
>> 3) Then  apply functions to account Dirichlet boundary condition:  
>>  VectorTools::interpolate_boundary_values(dof_handler, boundary_id, 
>> DirichletHandler, boundary_values);
>>
> ^^ that you would use during your setup
>  
>
>> MatrixTools::apply_boundary_values (boundary_values,
>> system_matrix,
>> solution,
>> system_rhs);
>>
>
> ^^ this is already what constraints.distribute_local_to_global() does, so 
> no need to do that.
>
> Regards
> Denis.
>  
>
>>
>>
>> Can you give some advice how to make it correct?
>>
>>
>> 
>>  
>> 
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>> On Tuesday, February 20, 2018 at 10:51:17 PM UTC+3, Denis Davydov wrote:
>>>
>>> Hi Ilya,
>>>
>>> I am not away of principally simple solution, the RHS simply contains 
>>> your bilinear form times interpolated values at constrained DoFs. 
>>> What you can do is to loop over a subset of cells (those which are at 
>>> the boundary), still assemble the local matrix and do 
>>> constraints.distribute_local_to_global() which takes global matrix and 
>>> RHS.
>>>
>>> Or you can use constrained_linear_operator 
>>> ,
>>>  
>>> see this 
>>>  forum 
>>> discussion.
>>>
>>> Please also study this 
>>>  
>>> thread, which covers similar question.
>>>
>>> Regards,
>>> Denis.
>>>
>>> On Tuesday, February 20, 2018 at 4:46:19 PM UTC+1, Bryukhanov Ilya wrote:

 Hi,

 I consider elastic problem and forces and displacements on boundaries 
 are changing with time. I don't want to compute FEM matrix
 on each step. How to correctly simulate the task? I know that there is 
 a simple solution, but my solution is incorrect :)
 I have two variables for the system_matrix and rhs_vector. I try to 
 recompute rhs_vector on each step  and 

[deal.II] Re: increment Dirichlet boundary condition

[Denis Davydov]

Hi Iliya,

there are numerous things that are wrong

On Wednesday, February 21, 2018 at 5:32:55 PM UTC+1, Bryukhanov Ilya wrote:
>
> Denis, thanks a lot for your answer and links that you provided. 
>
> However, my solution becomes in a way that the displacement values on that 
> boundary abnormally grow
> even without assigning new boundary condition in a time loop. I attach 
> figures on the zero step (correct) and on the first step (wrong).
>
> I do the following these things on each iteration: 
> 1) recompute rhs vector - system_rhs (which is zero for my problem),
> 2) compute matrix terms that are on the dirichlet boundary:
> for (auto cell = dof_handler.begin_active(); cell!=dof_handler.end(); 
> ++cell) {
> cell_matrix = 0; cell_rhs = 0; fe_values.reinit(cell);
> for (unsigned int face_number=0; 
> face_number if (cell->face(face_number)->at_boundary()) {
> // boundary index
> int id = cell->face(face_number)->boundary_id();
> if (boundary_displacement.find(id) != 
> boundary_displacement.end()) {
> //compute cell_matrix
>
> Then call function
> hanging_node_constraints.distribute_local_to_global(cell_matrix, cell_rhs, 
> local_dof_indices, system_matrix, system_rhs);
>

^^ that should be hanging_node_constraints + Dirichlet 

>
> 3) Then  apply functions to account Dirichlet boundary condition:  
>  VectorTools::interpolate_boundary_values(dof_handler, boundary_id, 
> DirichletHandler, boundary_values);
>
^^ that you would use during your setup
 

> MatrixTools::apply_boundary_values (boundary_values,
> system_matrix,
> solution,
> system_rhs);
>

^^ this is already what constraints.distribute_local_to_global() does, so 
no need to do that.

Regards
Denis.
 

>
>
> Can you give some advice how to make it correct?
>
>
> 
>  
> 
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
> On Tuesday, February 20, 2018 at 10:51:17 PM UTC+3, Denis Davydov wrote:
>>
>> Hi Ilya,
>>
>> I am not away of principally simple solution, the RHS simply contains 
>> your bilinear form times interpolated values at constrained DoFs. 
>> What you can do is to loop over a subset of cells (those which are at the 
>> boundary), still assemble the local matrix and do 
>> constraints.distribute_local_to_global() which takes global matrix and 
>> RHS.
>>
>> Or you can use constrained_linear_operator 
>> ,
>>  
>> see this 
>>  forum 
>> discussion.
>>
>> Please also study this 
>>  
>> thread, which covers similar question.
>>
>> Regards,
>> Denis.
>>
>> On Tuesday, February 20, 2018 at 4:46:19 PM UTC+1, Bryukhanov Ilya wrote:
>>>
>>> Hi,
>>>
>>> I consider elastic problem and forces and displacements on boundaries 
>>> are changing with time. I don't want to compute FEM matrix
>>> on each step. How to correctly simulate the task? I know that there is a 
>>> simple solution, but my solution is incorrect :)
>>> I have two variables for the system_matrix and rhs_vector. I try to 
>>> recompute rhs_vector on each step  and dirichlet condition 
>>> is accounted by calling "intepolate_boundary values + 
>>> apply_boundary_values" on each step with system_matrix and rhs_vector 
>>> variables inside.
>>> I think that I should use initial matrix and the matrix which is 
>>> modified by Dirichlet procedure and assign latter to the initial matrix on 
>>> each step?
>>> Am I corrent or there is some other solution?
>>>   
>>>
>>>
>>>
>>>
>>>
>>>

-- 
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 dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: increment Dirichlet boundary condition

[Bryukhanov Ilya]

Denis, thanks a lot for your answer and links that you provided. 

However, my solution becomes in a way that the displacement values on that 
boundary abnormally grow
even without assigning new boundary condition in a time loop. I attach 
figures on the zero step (correct) and on the first step (wrong).

I do the following these things on each iteration: 
1) recompute rhs vector - system_rhs (which is zero for my problem),
2) compute matrix terms that are on the dirichlet boundary:
for (auto cell = dof_handler.begin_active(); cell!=dof_handler.end(); 
++cell) {
cell_matrix = 0; cell_rhs = 0; fe_values.reinit(cell);
for (unsigned int face_number=0; 
face_numberface(face_number)->at_boundary()) {
// boundary index
int id = cell->face(face_number)->boundary_id();
if (boundary_displacement.find(id) != 
boundary_displacement.end()) {
//compute cell_matrix

Then call function
hanging_node_constraints.distribute_local_to_global(cell_matrix, cell_rhs, 
local_dof_indices, system_matrix, system_rhs);

3) Then  apply functions to account Dirichlet boundary condition:  
 VectorTools::interpolate_boundary_values(dof_handler, boundary_id, 
DirichletHandler, boundary_values);
MatrixTools::apply_boundary_values (boundary_values,
system_matrix,
solution,
system_rhs);



Can you give some advice how to make it correct?


 
















On Tuesday, February 20, 2018 at 10:51:17 PM UTC+3, Denis Davydov wrote:
>
> Hi Ilya,
>
> I am not away of principally simple solution, the RHS simply contains your 
> bilinear form times interpolated values at constrained DoFs. 
> What you can do is to loop over a subset of cells (those which are at the 
> boundary), still assemble the local matrix and do 
> constraints.distribute_local_to_global() which takes global matrix and RHS.
>
> Or you can use constrained_linear_operator 
> ,
>  
> see this  
> forum discussion.
>
> Please also study this 
>  thread, 
> which covers similar question.
>
> Regards,
> Denis.
>
> On Tuesday, February 20, 2018 at 4:46:19 PM UTC+1, Bryukhanov Ilya wrote:
>>
>> Hi,
>>
>> I consider elastic problem and forces and displacements on boundaries are 
>> changing with time. I don't want to compute FEM matrix
>> on each step. How to correctly simulate the task? I know that there is a 
>> simple solution, but my solution is incorrect :)
>> I have two variables for the system_matrix and rhs_vector. I try to 
>> recompute rhs_vector on each step  and dirichlet condition 
>> is accounted by calling "intepolate_boundary values + 
>> apply_boundary_values" on each step with system_matrix and rhs_vector 
>> variables inside.
>> I think that I should use initial matrix and the matrix which is modified 
>> by Dirichlet procedure and assign latter to the initial matrix on each step?
>> Am I corrent or there is some other solution?
>>   
>>
>>
>>
>>
>>
>>

-- 
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 dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

[deal.II] Re: increment Dirichlet boundary condition

[Denis Davydov]

Hi Ilya,

I am not away of principally simple solution, the RHS simply contains your 
bilinear form times interpolated values at constrained DoFs. 
What you can do is to loop over a subset of cells (those which are at the 
boundary), still assemble the local matrix and do 
constraints.distribute_local_to_global() which takes global matrix and RHS.

Or you can use constrained_linear_operator 
,
 
see this  
forum discussion.

Please also study this 
 thread, 
which covers similar question.

Regards,
Denis.

On Tuesday, February 20, 2018 at 4:46:19 PM UTC+1, Bryukhanov Ilya wrote:
>
> Hi,
>
> I consider elastic problem and forces and displacements on boundaries are 
> changing with time. I don't want to compute FEM matrix
> on each step. How to correctly simulate the task? I know that there is a 
> simple solution, but my solution is incorrect :)
> I have two variables for the system_matrix and rhs_vector. I try to 
> recompute rhs_vector on each step  and dirichlet condition 
> is accounted by calling "intepolate_boundary values + 
> apply_boundary_values" on each step with system_matrix and rhs_vector 
> variables inside.
> I think that I should use initial matrix and the matrix which is modified 
> by Dirichlet procedure and assign latter to the initial matrix on each step?
> Am I corrent or there is some other solution?
>   
>
>
>
>
>
>

-- 
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 dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.