Re: [deal.II] Re: imposing boundary conditions for a three-fields formulation

[Jean-Paul Pelteret]

Dear Alberto,

Now I understand your concern. The boolean flag that is passed as the final 
parameter to MatrixTools::apply_boundary_values 

 
determines whether column elimination is performed. You've set it to false, 
so this is not done and the RHS is unaltered. So switching it to true (the 
default) should preserve the symmetry of the system and make the necessary 
adjustments to the RHS vector.

J-P

On Thursday, April 20, 2017 at 2:19:50 AM UTC+2, Alberto Salvadori wrote:
>
> Dear Jean-Paul,
>
> thanks for your reply. 
>
> I see your point. It was my understanding (wrong, I guess) that in 
> imposing Dirichlet BC deal.II would modify the rows but also the columns, 
> adjusting the rhs accordingly. In this way, the symmetry of the system 
> matrix would have been preserved. I see now that this is not happening, in 
> all cases, is it correct? Then my question is about the solver, can we 
> still use CG as a solver if the system matrix looses symmetry? Does deal.II 
> make some further elaboration?
>
> Thank you again
>
> Alberto
>
>
>
>
>
> *Alberto Salvadori* Dipartimento di Ingegneria Civile, Architettura, 
> Territorio, Ambiente e di Matematica (DICATAM)
>  Universita` di Brescia, via Branze 43, 25123 Brescia
>  Italy
>  tel 030 3711239
>  fax 030 3711312
>
> e-mail: 
>  alberto.salvad...@unibs.it
> web-pages:
>  http://m4lab.unibs.it/faculty.html
>  http://dicata.ing.unibs.it/salvadori
>
> On Tue, Apr 18, 2017 at 11:10 AM, Jean-Paul Pelteret <> wrote:
>
>> Dear Alberto,
>>
>> The correctness of my answer depends on the numbering assigned to the 
>> degrees-of-freedom, which is something that you've not mentioned above. But 
>> from what I can discern the 9th row relates to the pressure test function, 
>> so K_{8,[0-7]} are the components from the coupling block K_{pu}, K_{88} is 
>> K_{pp} and K_{89} is K_{pJ}. Since these degree's of freedom are not 
>> constrained, this row should remain untouched after boundary values are 
>> applied. So, to me, nothing necessarily appears amiss on row 9 (as you've 
>> indicated).
>>
>> Do you agree with this?
>>
>> Kind regards,
>> Jean-Paul
>>
>> On Tuesday, April 18, 2017 at 12:32:44 PM UTC+2, Alberto Salvadori wrote:
>>>
>>> Hi, 
>>>
>>> sorry for asking your help again, that I really appreciated. 
>>> Today's issue is about imposing boundary conditions for a three-fields 
>>> formulation.
>>>
>>> In brief, I defined a 4 nodes square element, of type:
>>>
>>> fe(FE_Q(poly_degree), dim, // displacement
>>>
>>>FE_DGPMonomial(poly_degree - 1), 1, // pressure
>>>
>>>FE_DGPMonomial(poly_degree - 1), 1), // dilatation
>>>
>>> with  poly_degree = 1.
>>>
>>> While debugging on a single element, which has 10 dofs in 2D ( 4x2 
>>> displacements + 1 pressure + 1 dilatation), I aim at imposing Dirichlet 
>>> boundary conditions on displacements, whereas body forces are zero. 
>>>
>>> The cell matrix and the cell_rhs are the following.
>>>
>>> *Cell matrix = *
>>>
>>> *{  0.7037,   0.0278,  -0.2037,  -0.4722,  -0.1481,   0.4722,  -0.3519,  
>>> -0.0278,  -1.,   0.}, *
>>>
>>> *{  0.0278,   0.7037,   0.4722,  -0.1481,  -0.4722,  -0.2037,  -0.0278,  
>>> -0.3519,  -1.,   0.}, *
>>>
>>> *{ -0.2037,   0.4722,   0.7037,  -0.0278,  -0.3519,   0.0278,  -0.1481,  
>>> -0.4722,   1.,   0.}, *
>>>
>>> *{ -0.4722,  -0.1481,  -0.0278,   0.7037,   0.0278,  -0.3519,   0.4722,  
>>> -0.2037,  -1.,   0.}, *
>>>
>>> *{ -0.1481,  -0.4722,  -0.3519,   0.0278,   0.7037,  -0.0278,  -0.2037, 
>>>   0.4722,  -1.,   0.}, *
>>>
>>> *{  0.4722,  -0.2037,   0.0278,  -0.3519,  -0.0278,   0.7037,  -0.4722,  
>>> -0.1481,   1.,   0.}, *
>>>
>>> *{ -0.3519,  -0.0278,  -0.1481,   0.4722,  -0.2037,  -0.4722,   0.7037, 
>>>   0.0278,   1.,   0.}, *
>>>
>>> *{ -0.0278,  -0.3519,  -0.4722,  -0.2037,   0.4722,  -0.1481,   0.0278, 
>>>   0.7037,   1.,   0.}, *
>>>
>>> *{ -1.,  -1.,   1.,  -1.,  -1.,   1.,   1., 
>>>   1.,   0.,  -4.}, *
>>>
>>> *{  0.,   0.,   0.,   0.,   0.,   0.,   0., 
>>>   0.,  -4.,   6.6667}, *
>>>
>>> *Cell rhs =0.000e+00 0.000e+00 0.000e+00 0.000e+00 0.000e+00 0.000e+00 
>>> 0.000e+00 0.000e+00 0.000e+00 0.000e+00 *
>>>
>>> They seem to be correct. After calling
>>>
>>>   system_matrix.compress(VectorOperation::add);
>>>
>>>   system_rhs.compress(VectorOperation::add);
>>>
>>> the system matrix coincides with the cell matrix, the same for the rhs. 
>>> To impose BCs I followed the scheme of step-18, as follows:
>>>
>>>   std::map boundary_values;
>>>
>>>
>>>   ComponentMask displacements_mask = fe.component_mask (displacements);
>>>
>>>   VectorTools::interpolate_boundary_values (dof_handler,
>>>
>>> 0,
>>>
>>> 

Re: [deal.II] Re: imposing boundary conditions for a three-fields formulation

[Alberto Salvadori]

Dear Jean-Paul,

thanks for your reply.

I see your point. It was my understanding (wrong, I guess) that in imposing
Dirichlet BC deal.II would modify the rows but also the columns, adjusting
the rhs accordingly. In this way, the symmetry of the system matrix would
have been preserved. I see now that this is not happening, in all cases, is
it correct? Then my question is about the solver, can we still use CG as a
solver if the system matrix looses symmetry? Does deal.II make some further
elaboration?

Thank you again

Alberto





*Alberto Salvadori* Dipartimento di Ingegneria Civile, Architettura,
Territorio, Ambiente e di Matematica (DICATAM)
 Universita` di Brescia, via Branze 43, 25123 Brescia
 Italy
 tel 030 3711239
 fax 030 3711312

e-mail:
 alberto.salvad...@unibs.it
web-pages:
 http://m4lab.unibs.it/faculty.html
 http://dicata.ing.unibs.it/salvadori

On Tue, Apr 18, 2017 at 11:10 AM, Jean-Paul Pelteret 
wrote:

> Dear Alberto,
>
> The correctness of my answer depends on the numbering assigned to the
> degrees-of-freedom, which is something that you've not mentioned above. But
> from what I can discern the 9th row relates to the pressure test function,
> so K_{8,[0-7]} are the components from the coupling block K_{pu}, K_{88} is
> K_{pp} and K_{89} is K_{pJ}. Since these degree's of freedom are not
> constrained, this row should remain untouched after boundary values are
> applied. So, to me, nothing necessarily appears amiss on row 9 (as you've
> indicated).
>
> Do you agree with this?
>
> Kind regards,
> Jean-Paul
>
> On Tuesday, April 18, 2017 at 12:32:44 PM UTC+2, Alberto Salvadori wrote:
>>
>> Hi,
>>
>> sorry for asking your help again, that I really appreciated.
>> Today's issue is about imposing boundary conditions for a three-fields
>> formulation.
>>
>> In brief, I defined a 4 nodes square element, of type:
>>
>> fe(FE_Q(poly_degree), dim, // displacement
>>
>>FE_DGPMonomial(poly_degree - 1), 1, // pressure
>>
>>FE_DGPMonomial(poly_degree - 1), 1), // dilatation
>>
>> with  poly_degree = 1.
>>
>> While debugging on a single element, which has 10 dofs in 2D ( 4x2
>> displacements + 1 pressure + 1 dilatation), I aim at imposing Dirichlet
>> boundary conditions on displacements, whereas body forces are zero.
>>
>> The cell matrix and the cell_rhs are the following.
>>
>> *Cell matrix = *
>>
>> *{  0.7037,   0.0278,  -0.2037,  -0.4722,  -0.1481,   0.4722,  -0.3519,
>> -0.0278,  -1.,   0.}, *
>>
>> *{  0.0278,   0.7037,   0.4722,  -0.1481,  -0.4722,  -0.2037,  -0.0278,
>> -0.3519,  -1.,   0.}, *
>>
>> *{ -0.2037,   0.4722,   0.7037,  -0.0278,  -0.3519,   0.0278,  -0.1481,
>> -0.4722,   1.,   0.}, *
>>
>> *{ -0.4722,  -0.1481,  -0.0278,   0.7037,   0.0278,  -0.3519,   0.4722,
>> -0.2037,  -1.,   0.}, *
>>
>> *{ -0.1481,  -0.4722,  -0.3519,   0.0278,   0.7037,  -0.0278,  -0.2037,
>> 0.4722,  -1.,   0.}, *
>>
>> *{  0.4722,  -0.2037,   0.0278,  -0.3519,  -0.0278,   0.7037,  -0.4722,
>> -0.1481,   1.,   0.}, *
>>
>> *{ -0.3519,  -0.0278,  -0.1481,   0.4722,  -0.2037,  -0.4722,   0.7037,
>> 0.0278,   1.,   0.}, *
>>
>> *{ -0.0278,  -0.3519,  -0.4722,  -0.2037,   0.4722,  -0.1481,   0.0278,
>> 0.7037,   1.,   0.}, *
>>
>> *{ -1.,  -1.,   1.,  -1.,  -1.,   1.,   1.,
>> 1.,   0.,  -4.}, *
>>
>> *{  0.,   0.,   0.,   0.,   0.,   0.,   0.,
>> 0.,  -4.,   6.6667}, *
>>
>> *Cell rhs =0.000e+00 0.000e+00 0.000e+00 0.000e+00 0.000e+00 0.000e+00
>> 0.000e+00 0.000e+00 0.000e+00 0.000e+00 *
>>
>> They seem to be correct. After calling
>>
>>   system_matrix.compress(VectorOperation::add);
>>
>>   system_rhs.compress(VectorOperation::add);
>>
>> the system matrix coincides with the cell matrix, the same for the rhs.
>> To impose BCs I followed the scheme of step-18, as follows:
>>
>>   std::map boundary_values;
>>
>>
>>   ComponentMask displacements_mask = fe.component_mask (displacements);
>>
>>   VectorTools::interpolate_boundary_values (dof_handler,
>>
>> 0,
>>
>> IncrementalDirichletBoundary
>> Values( TIDM, NRit, GPDofs ),
>>
>> boundary_values,
>>
>> displacements_mask);
>>
>>
>>
>>   PETScWrappers::MPI::Vector tmp (locally_owned_dofs,mpi_communicator);
>>
>>   MatrixTools::apply_boundary_values (boundary_values,
>>
>>   system_matrix,
>>
>>   tmp,
>>
>>   system_rhs,
>>
>>   false);
>>
>>   incremental_displacement = tmp;
>>
>> where the class IncrementalDirichletBoundaryValues implements the
>> Dirichlet BC and GPDofs = 4.
>> The boundary values after VectorTools::interpolate_boundary_values 

[deal.II] Building UMFPACK with VS 2015

[Bill Greene]

I am continuing to try to use DEAL with Visual Studio.

I have had some success but ran into a problem when I tried to include
the bundled
umfpack via the Cmake directive
-D DEAL_II_WITH_UMFPACK=ON

Deal II seemed to build without problems but when I tried to build
example step-56,
I got unsatisfied externals: amd_malloc, amd_free, etc.

After some investigation it appears that those functions are being defined as:
nm -g amd_global.obj
 B ?amd_calloc@@3P6APEAX_K0@ZEA
0008 B ?amd_free@@3P6AXPEAX@ZEA
0018 B ?amd_malloc@@3P6APEAX_K@ZEA
 D ?amd_printf@@3P6AHPEBDZZEA
0010 B ?amd_realloc@@3P6APEAXPEAX_K@ZEA

but referenced as, for example:
$ nm -g umf_malloc.obj
 T ?umf_i_malloc@@YAPEAXH_K@Z
 U amd_malloc

I don't know if this is a fix or simply a work-around, but I added an
extern"C" to the definitions
of those functions in amd_global.cc as follows:

#ifdef __cplusplus
extern "C" {
#endif
/* standard ANSI-C: */
void *(*amd_malloc) (size_t) = malloc ;
void (*amd_free) (void *) = free ;
void *(*amd_realloc) (void *, size_t) = realloc ;
void *(*amd_calloc) (size_t, size_t) = calloc ;
#ifdef __cplusplus
}

Seems to have solved the problem for me.

Bill

-- 
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: QProjector::project_to_face and Gauss Lobatto

[Daniel Arndt]

vganpat,


I am trying to solve a problem that has both the volume integral and a 
> surface integral.
>
You would normally just create a quadrature rule for the volume and the 
surface integral separately.
 

> If I am using a higher order polynomial approximation and use Gauss 
> Lobatto quadrature, 
> the project_to_face() function seems to generate new points on the faces. 
> But, with GLL points, 
> there are already some points on the face that also contribute to the 
> volume integral. 
> What is not clear to me is, through this projection operation, I will lose 
> the unique global number 
> for these coincident points (those from Quadrature and then the 
> SubQuadrature). 
> If I were to record responses at the Quadrature points, with this 
> duplication of nodes, 
> how can I uniquely capture the response. Any help is appreciated. 
>
What exactly are you trying to achieve? In general, the quadrature points 
do not coincide with the support points of the used finite element
on a cell. Why are you using Gauss-Lobatto quadrature? Can you explain a 
bit more what you are doing?

Best,
Daniel

-- 
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.