Re: [deal.II] internal compiler error: Bus error

2017-02-01 Thread Matthias Maier
A bus error with gcc-4.8.1 on x86_64 is rare...

Is it reproducible? I.e. if you run the compilation with

  $ make VERBOSE=1 -j1

you will see the full command invocation that lead to the internal
compiler error. Does it always stop at the same .cc file ?

Further, can you ask the administrators to update the compiler?
gcc-4.9.4 is a well tested and solid choice. (the newer gcc-5* and
gcc-6* variants might lead to an ABI nightmare if you do not recompile
every C++ library).

Best,
Matthias



On Wed, Feb  1, 2017, at 20:05 CST, Lev Karatun  wrote:

> Hi Wolfgang,
>
> I was first compiling with -j8 then tried again with -j1 but got the same
> error.
>
> Best regards,
> Lev Karatun.
>
> 2017-02-01 20:41 GMT-05:00 Wolfgang Bangerth :
>
>> On 02/01/2017 06:23 PM, Lev Karatun wrote:
>>
>>>
>>> I'm trying to install the latest version of dealII on a cluster now, got
>>> another error unfortunately.
>>> The output is below, and I attached the logs too. I would appreciate any
>>> help with this.
>>>
>>
>> Are you compiling with -j4 or similar? Does the same happen if you compile
>> with -j1? The kind of internal compiler error you have here seems to happen
>> most often if you run out of memory.
>>
>> Best
>>  W.
>>
>> --
>> 
>> Wolfgang Bangerth  email: bange...@colostate.edu
>>www: http://www.math.colostate.edu/~bangerth/
>>
>> --
>> The deal.II project is located at http://www.dealii.org/
>> For mailing list/forum options, see https://groups.google.com/d/fo
>> rum/dealii?hl=en
>> --- You received this message because you are subscribed to a topic in the
>> Google Groups "deal.II User Group" group.
>> To unsubscribe from this topic, visit https://groups.google.com/d/to
>> pic/dealii/aZvol6kHWJI/unsubscribe.
>> To unsubscribe from this group and all its topics, send an email to
>> dealii+unsubscr...@googlegroups.com.
>> For more options, visit https://groups.google.com/d/optout.
>>

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


Re: [deal.II] internal compiler error: Bus error

2017-02-01 Thread Lev Karatun
Hi Wolfgang,

I was first compiling with -j8 then tried again with -j1 but got the same
error.

Best regards,
Lev Karatun.

2017-02-01 20:41 GMT-05:00 Wolfgang Bangerth :

> On 02/01/2017 06:23 PM, Lev Karatun wrote:
>
>>
>> I'm trying to install the latest version of dealII on a cluster now, got
>> another error unfortunately.
>> The output is below, and I attached the logs too. I would appreciate any
>> help with this.
>>
>
> Are you compiling with -j4 or similar? Does the same happen if you compile
> with -j1? The kind of internal compiler error you have here seems to happen
> most often if you run out of memory.
>
> Best
>  W.
>
> --
> 
> Wolfgang Bangerth  email: bange...@colostate.edu
>www: http://www.math.colostate.edu/~bangerth/
>
> --
> The deal.II project is located at http://www.dealii.org/
> For mailing list/forum options, see https://groups.google.com/d/fo
> rum/dealii?hl=en
> --- You received this message because you are subscribed to a topic in the
> Google Groups "deal.II User Group" group.
> To unsubscribe from this topic, visit https://groups.google.com/d/to
> pic/dealii/aZvol6kHWJI/unsubscribe.
> To unsubscribe from this group and all its topics, send an email to
> dealii+unsubscr...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.
>

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


Re: [deal.II] internal compiler error: Bus error

2017-02-01 Thread Wolfgang Bangerth

On 02/01/2017 06:23 PM, Lev Karatun wrote:


I'm trying to install the latest version of dealII on a cluster now, got
another error unfortunately.
The output is below, and I attached the logs too. I would appreciate any
help with this.


Are you compiling with -j4 or similar? Does the same happen if you 
compile with -j1? The kind of internal compiler error you have here 
seems to happen most often if you run out of memory.


Best
 W.

--

Wolfgang Bangerth  email: bange...@colostate.edu
   www: http://www.math.colostate.edu/~bangerth/

--
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] internal compiler error: Bus error

2017-02-01 Thread Lev Karatun
Hi everyone,

I'm trying to install the latest version of dealII on a cluster now, got 
another error unfortunately. 
The output is below, and I attached the logs too. I would appreciate any 
help with this. 

> In file included from 
>> /home/r/russ/lkaratun/aspect/dealii/include/deal.II/dofs/d 
>> 
>> of_handler.h:27:0,
>
>  from 
>> /home/r/russ/lkaratun/aspect/dealii/include/deal.II/dofs/d 
>> 
>> of_accessor.h:22,
>
>  from 
>> /home/r/russ/lkaratun/aspect/dealii/source/fe/fe_nedelec.c 
>> c:23:
>
> /home/r/russ/lkaratun/aspect/dealii/include/deal.II/dofs/block_info.h:233:34: 
>> in 
>> ternal compiler error: Bus error
>
>  types::global_dof_index BlockInfo::renumber (const unsigned int i) const
>
>   ^
>
> 0x8d7a9f crash_signal
>
> /scinet/gpc/src/gcc/gcc-4.8.1/gcc/toplev.c:332
>
> 0x5d4dc0 pop_deferring_access_checks()
>
> /scinet/gpc/src/gcc/gcc-4.8.1/gcc/cp/semantics.c:183
>
> 0x585848 cp_parser_parse_definitely
>
> /scinet/gpc/src/gcc/gcc-4.8.1/gcc/cp/parser.c:23820
>
> 0x59cc86 cp_parser_declarator
>
> /scinet/gpc/src/gcc/gcc-4.8.1/gcc/cp/parser.c:16227
>
> 0x5a28a0 cp_parser_init_declarator
>
> /scinet/gpc/src/gcc/gcc-4.8.1/gcc/cp/parser.c:15837
>
> 0x5a351e cp_parser_simple_declaration
>
> /scinet/gpc/src/gcc/gcc-4.8.1/gcc/cp/parser.c:10670
>
> 0x5a5020 cp_parser_block_declaration
>
> /scinet/gpc/src/gcc/gcc-4.8.1/gcc/cp/parser.c:10551
>
> 0x5acd6b cp_parser_declaration
>
> /scinet/gpc/src/gcc/gcc-4.8.1/gcc/cp/parser.c:10448
>
> 0x5abaad cp_parser_declaration_seq_opt
>
> /scinet/gpc/src/gcc/gcc-4.8.1/gcc/cp/parser.c:10334
>
> 0x5abcac cp_parser_namespace_body
>
> /scinet/gpc/src/gcc/gcc-4.8.1/gcc/cp/parser.c:15152
>
> 0x5abcac cp_parser_namespace_definition
>
> /scinet/gpc/src/gcc/gcc-4.8.1/gcc/cp/parser.c:15133
>
> 0x5accdd cp_parser_declaration
>
> /scinet/gpc/src/gcc/gcc-4.8.1/gcc/cp/parser.c:10436
>
> 0x5abaad cp_parser_declaration_seq_opt
>
> /scinet/gpc/src/gcc/gcc-4.8.1/gcc/cp/parser.c:10334
>
> 0x5ad2ba cp_parser_translation_unit
>
> /scinet/gpc/src/gcc/gcc-4.8.1/gcc/cp/parser.c:3813
>
> 0x5ad2ba c_parse_file()
>
> /scinet/gpc/src/gcc/gcc-4.8.1/gcc/cp/parser.c:28338
>
> 0x642a34 c_common_parse_file()
>
> /scinet/gpc/src/gcc/gcc-4.8.1/gcc/c-family/c-opts.c:1046
>
> Please submit a full bug report,
>
> with preprocessed source if appropriate.
>
> Please include the complete backtrace with any bug report.
>
> See  for instructions.
>
> make[2]: *** [source/fe/CMakeFiles/obj_fe_release.dir/fe_nedelec.cc.o] 
>> Error 1
>
> make[1]: *** [source/fe/CMakeFiles/obj_fe_release.dir/all] Error 2
>
> make: *** [all] Error 2
>
>
>

-- 
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 configuration:
#CMAKE_BUILD_TYPE:   Release
#BUILD_SHARED_LIBS:  ON
#CMAKE_INSTALL_PREFIX:   /home/r/russ/lkaratun/aspect/dealii
#CMAKE_SOURCE_DIR:   /home/r/russ/lkaratun/aspect/dealii
#(version 8.5.0-pre, shortrev f1284e9)
#CMAKE_BINARY_DIR:   /home/r/russ/lkaratun/aspect/dealii/build
#CMAKE_CXX_COMPILER: GNU 4.8.1 on platform Linux x86_64
#/scinet/gpc/compilers/gcc-4.8.1/bin/g++
#
#  Configured Features (DEAL_II_ALLOW_BUNDLED = ON, DEAL_II_ALLOW_AUTODETECTION 
= ON):
#  ( DEAL_II_WITH_64BIT_INDICES = OFF )
#  ( DEAL_II_WITH_ARPACK = OFF )
#DEAL_II_WITH_BOOST set up with bundled packages
#DEAL_II_WITH_BZIP2 set up with external dependencies
#  ( DEAL_II_WITH_CUDA = OFF )
#DEAL_II_WITH_CXX11 = ON
#  ( DEAL_II_WITH_CXX14 = OFF )
#  ( DEAL_II_WITH_GSL = OFF )
#  ( DEAL_II_WITH_HDF5 = OFF )
#  ( DEAL_II_WITH_LAPACK = OFF )
#  ( DEAL_II_WITH_METIS = OFF )
#DEAL_II_WITH_MPI set up with external dependencies
#DEAL_II_WITH_MUPARSER set up with bundled packages
#  ( DEAL_II_WITH_NETCDF = OFF )
#  ( DEAL_II_WITH_OPENCASCADE = OFF )
#DEAL_II_WITH_P4EST set up with external dependencies
#  ( DEAL_II_WITH_PETSC = OFF )
#  ( DEAL_II_WITH_SLEPC = OFF )
#DEAL_II_WITH_THREADS set up with bundled 

Re: [deal.II] using VectorTools::interpolate_boundary_values :Simple verification of Stokes system with zero velocity (using step-22)

2017-02-01 Thread Wolfgang Bangerth

On 02/01/2017 11:24 AM, Jane Lee wrote:


where the boundary values are just 1.0 so that u is a constant vector field.
this still gives me the bizarre 22 factor in my pressures solution
(which is just linear in z - the resulting pressure solution is atleast
linear), but the pressure solution aside, both the x and y components of
the velocity are 0.0 ON the boundaries and 1.0 everywhere else.


Are you calling constraints.distribute on the solution vector after 
solving your linear system?


Best
 W.

--

Wolfgang Bangerth  email: bange...@colostate.edu
   www: http://www.math.colostate.edu/~bangerth/

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


Re: [deal.II] using VectorTools::interpolate_boundary_values :Simple verification of Stokes system with zero velocity (using step-22)

2017-02-01 Thread Jane Lee
Hi Wolfgang, thanks for a swift reply. 

and apologies, I wasn't very clear. 

when I'm imposing dirichlet conditions everywhere, i don't have a condition 
for pressure (those were two separate test cases - sorry i was unclear). 
and i'm using:
VectorTools::interpolate_boundary_values (dof_handler,
  0,
  BoundaryValues(),
  constraints,
  
fe.component_mask(velocities));

where the boundary values are just 1.0 so that u is a constant vector field.
this still gives me the bizarre 22 factor in my pressures solution (which 
is just linear in z - the resulting pressure solution is atleast linear), 
but the pressure solution aside, both the x and y components of the 
velocity are 0.0 ON the boundaries and 1.0 everywhere else.

Yup sorry that was a mistake on my part - my real problem has total stress 
conditions so I naively just related it to pressure which isn't correct. 

would appreciate any help with why my interpolating boundary values aren't 
quite working.

Thanks again



On Wednesday, February 1, 2017 at 1:11:12 PM UTC, Wolfgang Bangerth wrote:
>
>
> > I'm trying to solve a stokes system with some funny boundary conditions 
> > so I've stripped the code to make sure I was getting the correct answers 
> > to a very simple system as follows: 
> > governing equations are: 
> > div grad u + grad p = f 
> > div u = 0 
> > i'm setting u=0 on the boundaries, giving u=0 (or close enough) 
> everywhere. 
> > setting f = -(0,1), i should get a linear profile with pressure 
> > increasing with depth at the same rate going downwards. 
>
> Yes. But if you impose Dirichlet boundary values for the velocity all 
> around the boundary, then the pressure is only determined up to a 
> constant, i.e., the pressure should be a linear function but the offset 
> is not determined by the equation. 
>
>
> > I'm clearly not understanding how to use 
> > VectorTools::interpolate_boundary_values correctly, and  I would 
> > appreciate some clarification on it. 
> > 
> > I'm imposing a zero pressure at the top of my rectangle with height of 
> > 10, which means I should get 10 at the bottom for pressure. 
>
> You can't do that. If you have Dirichlet boundary conditions for the 
> velocity, then you can't impose anything on the pressure. In fact, for 
> the incompressible Stokes equations, you can never impose any pressure 
> boundary conditions. 
>
> There is a description of boundary conditions for the Stokes equations 
> in the introduction of step-22. Take a look there to see what you can 
> and cannot impose for the Stokes equations. 
>
> Best 
>   W. 
>
> -- 
>  
> Wolfgang Bangerth  email: bang...@colostate.edu 
>  
> www: http://www.math.colostate.edu/~bangerth/ 
>

-- 
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: Source of a constant residual error when a triangulation has hanging nodes

2017-02-01 Thread amir . pasha . beh
Dear Daniel

The code solve a solid mechanic model with a plastic material model. The 
WorkStream is used for matrix and residual assembly (like step-44). So I 
just add some new classes  (PerTaskData and ScratchData) and the most parts 
of the code is not changed. The code works similar to the its serial 
version without any problem when the local adaptive refinement is not 
active. During a load increment, let the solution to time step Tn be given. 
The full problem is solved. First, I refine the mesh based on the new 
solution and interpolate the old solution (at Tn ) onto the new mesh. Then 
I solve for the solution at Tn+1 again, but on the refined mesh. This 
process is tested in the previous version of the code. However, in the new 
code after the first refinement, the solution for Tn+ has a constant 
residual error during Newton iteration. 

Thank you
P.

On Wednesday, February 1, 2017 at 8:24:51 PM UTC+3:30, Daniel Arndt wrote:
>
> Pasha,
>
> I have a code that works fine. Now, I am working to change the code 
>> structure to use “workstream”. The new code works fine without adaptive 
>> refinement. However, when the mesh is refined and the triangulation has 
>> hanging node, a constant residual error remains during Newton iteration. I 
>> cant find source of the problem
>>
> You are providing way too few information to expect useful answers. What 
> are you doing in your code?
> For which part in your code do you want to change to WorkStream? 
> What exactly did you change? How is this related to whether you are doing 
> adaptive refinement or not?
>
> 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.


[deal.II] Re: Source of a constant residual error when a triangulation has hanging nodes

2017-02-01 Thread Daniel Arndt
Pasha,

I have a code that works fine. Now, I am working to change the code 
> structure to use “workstream”. The new code works fine without adaptive 
> refinement. However, when the mesh is refined and the triangulation has 
> hanging node, a constant residual error remains during Newton iteration. I 
> cant find source of the problem
>
You are providing way too few information to expect useful answers. What 
are you doing in your code?
For which part in your code do you want to change to WorkStream? 
What exactly did you change? How is this related to whether you are doing 
adaptive refinement or not?

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.


[deal.II] Re: Shape function derivatives

2017-02-01 Thread 'Seyed Ali Mohseni' via deal.II User Group
@Prof. Bangerth: Thank you very much :) Finally, I understand how the 
stiffness matrix in step-18 is constructed. Your shape function derivatives 
are basically already the B-operator. 

Kind regards,
S. A. Mohseni

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


Re: [deal.II] Shape function derivatives

2017-02-01 Thread Wolfgang Bangerth

On 02/01/2017 07:02 AM, 'Seyed Ali Mohseni' via deal.II User Group wrote:

My approach is the following:

|
template
FullMatrix get_b_operator (const FEValues _values, const
unsigned int dofs_per_cell, const unsigned int q)
{
FullMatrix tmp(dim, GeometryInfo::vertices_per_cell);

std::vector > invJ =
fe_values.get_inverse_jacobians();

for (unsigned int i = 0; i < dofs_per_cell; i = i + dim)
{
const unsigned int index = i / dim;

// COMPUTE: B-operator (Remark: This version has to be extended for 3D!)
tmp[0][index] = invJ[q][0][0] * fe_values.shape_grad(i, q)[0] +
invJ[q][0][1] * fe_values.shape_grad(i, q)[1];
tmp[1][index] = invJ[q][1][0] * fe_values.shape_grad(i, q)[0] +
invJ[q][1][1] * fe_values.shape_grad(i, q)[1];
}

return tmp;
}
|


fe_values.shape_grad() returns the gradient *on the actual cell*, not on 
the reference cell. Consequently, you do not have to multiply it by invJ 
again. This is the source of your factor of two.


Best
 W.

--

Wolfgang Bangerth  email: bange...@colostate.edu
   www: http://www.math.colostate.edu/~bangerth/

--
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] Shape function derivatives

2017-02-01 Thread 'Seyed Ali Mohseni' via deal.II User Group
Hi,

I found some issue when constructing the FE B-operator for multiple 
elements. In a single element it works perfectly.
My approach is the following:

template
FullMatrix get_b_operator (const FEValues _values, const 
unsigned int dofs_per_cell, const unsigned int q)
{
FullMatrix tmp(dim, GeometryInfo::vertices_per_cell);

std::vector > invJ = 
fe_values.get_inverse_jacobians();

for (unsigned int i = 0; i < dofs_per_cell; i = i + dim)
{
const unsigned int index = i / dim;

// COMPUTE: B-operator (Remark: This version has to be extended for 3D!)
tmp[0][index] = invJ[q][0][0] * fe_values.shape_grad(i, q)[0] + 
invJ[q][0][1] * fe_values.shape_grad(i, q)[1];
tmp[1][index] = invJ[q][1][0] * fe_values.shape_grad(i, q)[0] + 
invJ[q][1][1] * fe_values.shape_grad(i, q)[1];
}

return tmp;
}

Hence, I am computing B = inv(J) * transpose(dN). For a single element I 
receive the following correct results for dN:

*C++ CODE*


SHAPE FUNCTION DERIVATIVES (GAUSS POINT 1)

-0.788675  -0.788675  
0.788675  -0.211325  
0.211325  0.211325  
-0.211325  0.788675  

SHAPE FUNCTION DERIVATIVES (GAUSS POINT 2)

-0.788675  -0.211325  
0.788675  -0.788675  
0.211325  0.788675  
-0.211325  0.211325  

SHAPE FUNCTION DERIVATIVES (GAUSS POINT 3)

-0.211325  -0.788675  
0.211325  -0.211325  
0.788675  0.211325  
-0.788675  0.788675  

SHAPE FUNCTION DERIVATIVES (GAUSS POINT 4)
-0.211325  -0.211325  
0.211325  -0.788675  
0.788675  0.788675  
-0.788675  0.211325  



*DEAL.II*


SHAPE FUNCTION DERIVATIVES

  -0.788675 -0.788675   -0.788675 -0.788675   0.788675 -0.211325   0.788675 
-0.211325   -0.211325 0.788675   -0.211325 0.788675   0.211325 0.211325   
0.211325 0.211325 
  -0.788675 -0.211325   -0.788675 -0.211325   0.788675 -0.788675   0.788675 
-0.788675   -0.211325 0.211325   -0.211325 0.211325   0.211325 0.788675   
0.211325 0.788675 
  -0.211325 -0.788675   -0.211325 -0.788675   0.211325 -0.211325   0.211325 
-0.211325   -0.788675 0.788675   -0.788675 0.788675   0.788675 0.211325   
0.788675 0.211325 
  -0.211325 -0.211325   -0.211325 -0.211325   0.211325 -0.788675   0.211325 
-0.788675   -0.788675 0.211325   -0.788675 0.211325   0.788675 0.788675   
0.788675 0.788675 


Due to the ordering of the element connectivity (numbering of nodes are 
different in deal.II) we have to swap column 3 and 4 for comparison.

Now for 4 elements there is something strange when I compute my shape 
function derivatives within deal.II:

*DEAL.II*


SHAPE FUNCTION DERIVATIVES

  -1.577350 -1.577350   -1.577350 -1.577350   1.577350 -0.422650   1.577350 
-0.422650   -0.422650 1.577350   -0.422650 1.577350   0.422650 0.422650   
0.422650 0.422650 
  -1.577350 -0.422650   -1.577350 -0.422650   1.577350 -1.577350   1.577350 
-1.577350   -0.422650 0.422650   -0.422650 0.422650   0.422650 1.577350   
0.422650 1.577350 
  -0.422650 -1.577350   -0.422650 -1.577350   0.422650 -0.422650   0.422650 
-0.422650   -1.577350 1.577350   -1.577350 1.577350   1.577350 0.422650   
1.577350 0.422650 
  -0.422650 -0.422650   -0.422650 -0.422650   0.422650 -1.577350   0.422650 
-1.577350   -1.577350 0.422650   -1.577350 0.422650   1.577350 1.577350   
1.577350 1.577350 



 The results for one element change now by a factor of 2 and with 
increasing element size it always differs by another factor of 2 again.

I was wondering, if my approach to compute my B-operator was wrong. Doesn't 
fe_values.shape_grad() give me the shape function derivatives directly? Or 
do I have to use shape_grad_component? 

Additionally, I checked the symmetric_gradient from step-18 and it contains 
actually the correct values I need, but somehow in a strange order:


SYMMETRIC GRADIENT (GAUSS POINT 1)

0.00 0.211325 
0.211325 0.422650

SYMMETRIC GRADIENT (GAUSS POINT 2)

0.00 0.211325 
0.211325 1.577350

SYMMETRIC GRADIENT (GAUSS POINT 3)

0.00 0.788675 
0.788675 0.422650 

SYMMETRIC GRADIENT (GAUSS POINT 4)

0.00 0.788675 
0.788675 1.577350




Kind regards,
S. A. Mohseni


-- 
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] Source of a constant residual error when a triangulation has hanging nodes

2017-02-01 Thread amir . pasha . beh
 

Dear all


I have a code that works fine. Now, I am working to change the code 
structure to use “workstream”. The new code works fine without adaptive 
refinement. However, when the mesh is refined and the triangulation has 
hanging node, a constant residual error remains during Newton iteration. I 
cant find source of the problem. 

There are two “constraints matrix”. On of them for “Newton update” and one 
for “hanging nodes”. The constraints matrix are merged and is used to 
determine error as:


BlockVector<*double*> error_ud(dofs_per_block);

*for* (*unsigned* *int* i = 0; i < dof_handler_ref.n_dofs(); ++i)

*if* (!constraints_update.*is_constrained*(i))

error_ud(i) = system_rhs(i);


error_update.norm = error_ud.l2_norm();



Regards

Pasha

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


Re: [deal.II] using VectorTools::interpolate_boundary_values :Simple verification of Stokes system with zero velocity (using step-22)

2017-02-01 Thread Wolfgang Bangerth



I'm trying to solve a stokes system with some funny boundary conditions
so I've stripped the code to make sure I was getting the correct answers
to a very simple system as follows:
governing equations are:
div grad u + grad p = f
div u = 0
i'm setting u=0 on the boundaries, giving u=0 (or close enough) everywhere.
setting f = -(0,1), i should get a linear profile with pressure
increasing with depth at the same rate going downwards.


Yes. But if you impose Dirichlet boundary values for the velocity all 
around the boundary, then the pressure is only determined up to a 
constant, i.e., the pressure should be a linear function but the offset 
is not determined by the equation.




I'm clearly not understanding how to use
VectorTools::interpolate_boundary_values correctly, and  I would
appreciate some clarification on it.

I'm imposing a zero pressure at the top of my rectangle with height of
10, which means I should get 10 at the bottom for pressure.


You can't do that. If you have Dirichlet boundary conditions for the 
velocity, then you can't impose anything on the pressure. In fact, for 
the incompressible Stokes equations, you can never impose any pressure 
boundary conditions.


There is a description of boundary conditions for the Stokes equations 
in the introduction of step-22. Take a look there to see what you can 
and cannot impose for the Stokes equations.


Best
 W.

--

Wolfgang Bangerth  email: bange...@colostate.edu
   www: http://www.math.colostate.edu/~bangerth/

--
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] using VectorTools::interpolate_boundary_values :Simple verification of Stokes system with zero velocity (using step-22)

2017-02-01 Thread Jane Lee
Hi all,

This may seem like a very obvious question/mistake on my part so please do 
excuse me if I've overseen something very simple. 

I'm trying to solve a stokes system with some funny boundary conditions so 
I've stripped the code to make sure I was getting the correct answers to a 
very simple system as follows:
governing equations are:
div grad u + grad p = f
div u = 0
i'm setting u=0 on the boundaries, giving u=0 (or close enough) everywhere.
setting f = -(0,1), i should get a linear profile with pressure increasing 
with depth at the same rate going downwards. 

I'm clearly not understanding how to use 
VectorTools::interpolate_boundary_values correctly, and  I would appreciate 
some clarification on it. 

I'm imposing a zero pressure at the top of my rectangle with height of 10, 
which means I should get 10 at the bottom for pressure. 
I'm using VectorTools::interpolate_boundary_values with component mask to 
do this, like in step-22 but I seem to get an odd result, where there is a 
factor of 22 that is multiplied for the pressure. 
This way just generally doesn't seem to work, eg. if i set the boundary 
value as 10.0 using the same way, I don't get 10.0 for pressure at the top, 
I still get 0.0 for the top, but at the bottom, instead of 220 which i get 
with 0.0 condition at the top, i get 230.
Similarly, when i try (or think i am) setting u=1 on the boundaries (and 
therefore u=1 everywhere), the output is showing 1 everywhere EXCEPT on the 
boundaries, which seems very odd

Could anyone possibly tell me what I'm doing wrong to with interpolating 
the boundary values? seems like such a simple issue but i can't seem to 
figure it out

code is attached - boundary id 0 are the sides, 1 the top and 2 the bottom. 
NB i am just using a direct solver for now - the rest is taken mostly from 
step-22.

Many thanks

-- 
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.
/* -
 *
 * Copyright (C) 2008 - 2015 by the deal.II authors
 *
 * This file is part of the deal.II library.
 *
 * The deal.II library is free software; you can use it, redistribute
 * it, and/or modify it under the terms of the GNU Lesser General
 * Public License as published by the Free Software Foundation; either
 * version 2.1 of the License, or (at your option) any later version.
 * The full text of the license can be found in the file LICENSE at
 * the top level of the deal.II distribution.
 *
 * -
 
 *
 * Author: Wolfgang Bangerth, Texas A University, 2008
 */



#include 
#include 
#include 
#include 

#include 
#include 
#include 
#include 
#include 
#include 

#include 
#include 
#include 
#include 
#include 
#include 
#include 

#include 
#include 
#include 
#include 

#include 
#include 
#include 

#include 
#include 
#include 
#include 

#include 

#include 

#include 
#include 
#include 

namespace Step22
{
using namespace dealii;

namespace data
{
const double rho_B = 1.0;
const double eta = 1.0;
const double top = 10.0;
const double bottom = 0.0;
const int dimension = 2;

}




template 
class StokesProblem
{
public:
StokesProblem (const unsigned int degree);
void run ();

private:
void setup_dofs ();
void assemble_system ();
void solve ();
void output_results (const unsigned int refinement_cycle) const;
void refine_mesh ();

const unsigned int   degree;

Triangulation   triangulation;
FESystemfe;
DoFHandler  dof_handler;

ConstraintMatrix constraints;

BlockSparsityPattern  sparsity_pattern;
BlockSparseMatrix system_matrix;

BlockVector solution;
BlockVector system_rhs;


};



template 
class BoundaryValues : public Function
{
public:
  BoundaryValues () : Function(dim+1) {}
  virtual double value (const Point   ,
const unsigned int  component = 0) const;
  virtual void vector_value (const Point ,
 Vector   ) const;
};
template 
double
BoundaryValues::value (const Point  ,
const unsigned int component) const
{

  return 1.0;
}
template 
void
BoundaryValues::vector_value (const Point ,
   Vector   ) const
{