Re: [deal.II] Diffusion equation with DG fails with distorted mesh
On 4/25/19 2:50 PM, Gary Uppal wrote: > > I am trying to solve the diffusion equation with Discontinuous Galerkin > elements. The solution looks good with a regular structured mesh, but if I > distort the mesh, the solution blows up and does not converge. Is there an > obvious reason this would happen? I later need a mesh with holes in it, so I > cannot always use the structured mesh. > > I use the interior penalty method and get the diffusion matrix using > MeshWorker as in Step-39. I compute the mass matrix and solve the diffusion > equation with an implicit backward Euler method and am using periodic > boundary > conditions. Snapshots of the solution with structured and distorted meshes > are > shown below. Any help is appreciated! Gary -- the wrong solutions happen on cells that are very nearly degenerate (i.e., are almost triangular). The finite element theory says that the error between the exact and the numerical solution is bounded by a constant times some power of h, where the constant depends on the minimal and maximal angles at the vertices of the cells. Theory then also predicts that this constant goes to infinity if the maximal angle at one vertex comes close to 180 degrees -- which is exactly what is happening in your case. So choose a mesh that is less distorted and you should be fine. 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] Re: Integrating a deal.II-specific function in a NOX-Interface for MPI-threads
On 4/25/19 10:13 AM, 'Maxi Miller' via deal.II User Group wrote: > After running some tests, I ended up reducing the problem to the transfer to > and from the Epetra-vectors. I have to write an interface to the model > (according to the instructions), and as shown in the code above. There I have > Epetra-Vectors created by my interface, with a size of > locally_relevant_dofs.size(), and TrilinosWrappers::MPI::Vectors. Transfer > from the Epetra-Vectors to the MPI::Vectors works fine, but the way back does > not work (The MPI::Vectors are larger than the Epetra_Vectors). Are there any > hints in how I still could fit the data from the MPI::Vector into the > Epetra_Vector? Or should I rather ask this on the Trilinos mailing list? Good question. I think it would probably be very useful to have small testcase others could look at. The program you have attached is still 1,300 lines, which is far too much for anyone to understand. But since you have an idea of what the problem is, do you think you could come up with a small testcase that illustrates exactly the issue you have? It doesn't have to do anything useful at all -- in your case, I think all that's necessary is to create a vector, convert it as you describe above, and then output the sizes to show that they are not the same. Run this on two processors, and you should have something that could be done in 50 or 100 lines, and I think that might be small enough for someone who doesn't know your code to understand what is necessary. I've built these sorts of testcases either from scratch in the past, or by just keep removing things from a program that (i) are run after the time the problem happens, and (ii) then removing everything that is not necessary. 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] Setting an initial condition by fixing the DOF values themselves
Dear Wolfgang, this was easier than I thought. I followed your suggestion and used VectorTools::interpolate. Everything works as expected. Thank you very much Final solution: // initialConditionParameters.uvw is the parsed function // newton_update is a local vector const FEValuesExtractors::Vector velocities (0); const FEValuesExtractors::Scalar pressure (dim); const MappingQ mapping (degreeVelocity_,femParameters.qmapping_all); VectorTools::interpolate(mapping, dof_handler, initialConditionParameters.uvw, newton_update, fe.component_mask(velocities)); Thanks again! On Thursday, 25 April 2019 09:32:09 UTC-4, Wolfgang Bangerth wrote: > > On 4/25/19 6:58 AM, Bruno Blais wrote: > > I would like to be able to set-up an initial condition by fixing the > value of > > the DOF themselves by using their x,y(,z) position with a ParsedFuction. > > > > Generally, I set the initial condition by using an L2 projection of the > > function (for complex function and higher order element I know this is > more > > appropriate). > > However, for some stuff I would like to be able to fix the value of the > DOF > > directly. > > My issue is, through the doxygen, I have not found a way to loop through > the > > DOF and to extract their position, even more when I have high order > element or > > when I have a multiplicity (say velocity-vector + pressure). > > > > What I would be looking into doing would be something similar to : > > > > loop over local DOF; > >get DOF position > >get DOF component (u[0], u[1], u[2], p) > >calculate parsed function using the point > >set the value of the DOF > > end > > > > Is it something that is possible? This is a weird way of iterating since > I > > know we always loop over the cells. Maybe there is a way to achieve this > by > > looping over the cells and then looping over the DOF of the cells? This > would > > be slightly redundant, but I do not care since I am doing this only > once. > > You can get the locations of DoFs using > DoFTools::map_dofs_to_support_points(). > > But maybe the easier way is to use > VectorTools::interpolate_boundary_values() > or project() and just pass it a function object that describes the > function > you want to set in terms of (x,y,z). This could be your ParsedFunction. > interpolate_b_v() will then call your parsed function with a bunch of > points > that happen to correspond to the DoF locations, but it makes the > connection > internally and you no longer need to worry about which DoF sits where. > > 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.
Re: [deal.II] Setting an initial condition by fixing the DOF values themselves
Dear Wolfgang, This was easier than I thought. I followed your suggestion and used VectorTools::interpolate Final solution to set the values of the DOF using a parsed function: // initialConditionParameters.uvw is my parsed function // newtonUpdate is my local vector that I use to alter the DOF values const FEValuesExtractors::Vector velocities (0); const FEValuesExtractors::Scalar pressure (dim); const MappingQ mapping (degreeVelocity_,femParameters.qmapping_all); VectorTools::interpolate(mapping, dof_handler, initialConditionParameters.uvw, newton_update, fe.component_mask(velocities)); Thank you very much for the help! On Thursday, 25 April 2019 09:32:09 UTC-4, Wolfgang Bangerth wrote: > > On 4/25/19 6:58 AM, Bruno Blais wrote: > > I would like to be able to set-up an initial condition by fixing the > value of > > the DOF themselves by using their x,y(,z) position with a ParsedFuction. > > > > Generally, I set the initial condition by using an L2 projection of the > > function (for complex function and higher order element I know this is > more > > appropriate). > > However, for some stuff I would like to be able to fix the value of the > DOF > > directly. > > My issue is, through the doxygen, I have not found a way to loop through > the > > DOF and to extract their position, even more when I have high order > element or > > when I have a multiplicity (say velocity-vector + pressure). > > > > What I would be looking into doing would be something similar to : > > > > loop over local DOF; > >get DOF position > >get DOF component (u[0], u[1], u[2], p) > >calculate parsed function using the point > >set the value of the DOF > > end > > > > Is it something that is possible? This is a weird way of iterating since > I > > know we always loop over the cells. Maybe there is a way to achieve this > by > > looping over the cells and then looping over the DOF of the cells? This > would > > be slightly redundant, but I do not care since I am doing this only > once. > > You can get the locations of DoFs using > DoFTools::map_dofs_to_support_points(). > > But maybe the easier way is to use > VectorTools::interpolate_boundary_values() > or project() and just pass it a function object that describes the > function > you want to set in terms of (x,y,z). This could be your ParsedFunction. > interpolate_b_v() will then call your parsed function with a bunch of > points > that happen to correspond to the DoF locations, but it makes the > connection > internally and you no longer need to worry about which DoF sits where. > > 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.
Re: [deal.II] Imposing the Neumann BC using right_hand_side() function
On 4/25/19 9:46 AM, Muhammad Mashhood wrote: > > As unlike the step-8 here in step-26 the rhs is being formed directly > globally as "system_rhs" using "create_right_hand_side()" function so > can I use the "create_boundary_right_hand_side()" function using > following lines of code? It is definitely possible to *add* to the vector you get from create_boundary_rhs(). But in this call... > / > / > / std::map::iterator > it=boundary_count.begin(); > > VectorTools::create_boundary_right_hand_side(dof_handler, > > QGauss(fe.degree+1), > rhs_function, > tmp,/ > / it->first);/ > where the "rhs_function" provides the values of /g0/ and /g1 /at > boundary node points and "it->first" is the boundary id where I want to > apply the BC. ...you are asked to pass in a Quadrature, not a Quadrature. Consequently, the code does not compile. This is actually stated quite explicitly in the second part of the error message: dealii::VectorTools::create_boundary_right_hand_side(const dealii::DoFHandler&, const dealii::Quadrature<(dim - 1)>&, const dealii::Function&, VectorType&, const std::set&) [with int dim = 1; int spacedim = 1; VectorType = dealii::Vector; typename VectorType::value_type = double] void create_boundary_right_hand_side ^ /usr/local/include/deal.II/numerics/vector_tools.h:2105:8: note: no known conversion for argument 2 from ‘dealii::QGauss<1>’ to ‘const dealii::Quadrature<0>&’ Best WB -- 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] Re: Integrating a deal.II-specific function in a NOX-Interface for MPI-threads
After running some tests, I ended up reducing the problem to the transfer
to and from the Epetra-vectors. I have to write an interface to the model
(according to the instructions), and as shown in the code above. There I
have Epetra-Vectors created by my interface, with a size of
locally_relevant_dofs.size(), and TrilinosWrappers::MPI::Vectors. Transfer
from the Epetra-Vectors to the MPI::Vectors works fine, but the way back
does not work (The MPI::Vectors are larger than the Epetra_Vectors). Are
there any hints in how I still could fit the data from the MPI::Vector into
the Epetra_Vector? Or should I rather ask this on the Trilinos mailing list?
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) 1999 - 2018 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.md at
* the top level directory of deal.II.
*
* -
*
* Author: Wolfgang Bangerth, University of Heidelberg, 1999
*/
// @sect3{Include files}
//Nox include files
#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
#include
#include
#include
#include
#include
#include
#include
// NOX Objects
#include "NOX.H"
#include "NOX_Epetra.H"
// Trilinos Objects
#ifdef HAVE_MPI
#include "Epetra_MpiComm.h"
#else
#include "Epetra_SerialComm.h"
#endif
#include "Epetra_Vector.h"
#include "Epetra_Operator.h"
#include "Epetra_RowMatrix.h"
#include "NOX_Epetra_Interface_Required.H" // base class
#include "NOX_Epetra_Interface_Jacobian.H" // base class
#include "NOX_Epetra_Interface_Preconditioner.H" // base class
#include "Epetra_CrsMatrix.h"
#include "Epetra_Map.h"
#include "Epetra_LinearProblem.h"
#include "AztecOO.h"
#include "Teuchos_StandardCatchMacros.hpp"
// User's application specific files
//#include "problem_interface.h" // Interface file to NOX
//#include "finite_element_problem.h"
// Required for reading and writing parameter lists from xml format
// Configure Trilinos with --enable-teuchos-extended
#ifdef HAVE_TEUCHOS_EXTENDED
#include "Teuchos_XMLParameterListHelpers.hpp"
#include
#endif
// This is C++ ...
#include
#include
using namespace dealii;
template
class BoundaryValues : public dealii::Function
{
public:
BoundaryValues(const size_t n_components)
: dealii::Function(1), n_components(n_components)
{}
virtual double value (const dealii::Point &p,
const unsigned int component = 0) const;
virtual void vector_value(const dealii::Point &p, dealii::Vector &value) const;
private:
const size_t n_components;
};
template
double BoundaryValues::value (const dealii::Point &p,
const unsigned int) const
{
const double x = p[0];
const double y = p[1];
return sin(M_PI * x) * sin(M_PI * y);
}
template
void BoundaryValues::vector_value(const dealii::Point &p, dealii::Vector &value) const
{
for(size_t i = 0; i < value.size(); ++i)
value[i] = BoundaryValues::value(p, i);
}
template
class Solution : public Function
{
public:
Solution() : Function(1)
{
}
virtual double value(const Point &p, const unsigned int component) const override;
virtual Tensor<1, dim> gradient(const Point &p, const unsigned int component) const override;
private:
};
template
double Solution::value(const Point &p, const unsigned int) const
{
const double x = p[0];
const double y = p[1];
return sin(M_PI * x) * sin(M_PI * y);
}
template
Tensor<1, dim> Solution::gradient(const Point &p, const unsigned int) const
{
Tensor<1, dim> return_value;
AssertThrow(dim == 2, ExcNotImplemented());
const double x = p[0];
const double y = p[1];
return_value[0] = M_PI * cos(M_PI * x) * sin(M_PI * y);
return_value[1] = M_PI * cos(M_PI * y) * sin(M_PI * x);
return return_value;
}
template
class ProblemIn
Re: [deal.II] Imposing the Neumann BC using right_hand_side() function
In the weak form after Neuman BC implementation for equation in step-26 i got following additional terms contributing in rhs: *for 1D domain x = [0, 1]* *system_rhs + k_n * [ g1 * phi_(x1) - g0 * phi_(x0) ], where g0 and g1 are the flux values at boundary nodes and phi_(x0), phi_(x1) are shape functions at these boundary nodes on both ends of a bar.* As unlike the step-8 here in step-26 the rhs is being formed directly globally as "system_rhs" using "create_right_hand_side()" function so can I use the "create_boundary_right_hand_side()" function using following lines of code? *std::map::iterator it=boundary_count.begin(); VectorTools::create_boundary_right_hand_side(dof_handler, QGauss(fe.degree+1), rhs_function, tmp,* * it->first);* where the "rhs_function" provides the values of *g0* and *g1 *at boundary node points and "it->first" is the boundary id where I want to apply the BC. So far I have tried it from my side but giving the following error: */home/muhammad/software/dealii-8.5.0/examples/step-26/step-26.cc:602:53: error: no matching function for call to ‘create_boundary_right_hand_side(dealii::DoFHandler<1>&, dealii::QGauss<1>, Step26::RightHandSide<1>&, dealii::Vector&)’ VectorTools::create_boundary_right_hand_side(dof_handler, ^In file included from /home/muhammad/software/dealii-8.5.0/examples/step-26/step-26.cc:46:0:/usr/local/include/deal.II/numerics/vector_tools.h:2090:8: note: candidate: template void dealii::VectorTools::create_boundary_right_hand_side(const dealii::Mapping&, const dealii::DoFHandler&, const dealii::Quadrature<(dim - 1)>&, const dealii::Function&, VectorType&, const std::set&) void create_boundary_right_hand_side (const Mapping&mapping, ^/usr/local/include/deal.II/numerics/vector_tools.h:2090:8: note: template argument deduction/substitution failed:/home/muhammad/software/dealii-8.5.0/examples/step-26/step-26.cc:602:53: note: ‘dealii::DoFHandler<1>’ is not derived from ‘const dealii::Mapping’ VectorTools::create_boundary_right_hand_side(dof_handler, ^In file included from /home/muhammad/software/dealii-8.5.0/examples/step-26/step-26.cc:46:0:/usr/local/include/deal.II/numerics/vector_tools.h:2105:8: note: candidate: void dealii::VectorTools::create_boundary_right_hand_side(const dealii::DoFHandler&, const dealii::Quadrature<(dim - 1)>&, const dealii::Function&, VectorType&, const std::set&) [with int dim = 1; int spacedim = 1; VectorType = dealii::Vector; typename VectorType::value_type = double] void create_boundary_right_hand_side ^/usr/local/include/deal.II/numerics/vector_tools.h:2105:8: note: no known conversion for argument 2 from ‘dealii::QGauss<1>’ to ‘const dealii::Quadrature<0>&’/usr/local/include/deal.II/numerics/vector_tools.h:2119:8: note: candidate: template void dealii::VectorTools::create_boundary_right_hand_side(const dealii::hp::MappingCollection&, const dealii::hp::DoFHandler&, const dealii::hp::QCollection<(dim - 1)>&, const dealii::Function&, VectorType&, const std::set&) void create_boundary_right_hand_side ^/usr/local/include/deal.II/numerics/vector_tools.h:2119:8: note: template argument deduction/substitution failed:/home/muhammad/software/dealii-8.5.0/examples/step-26/step-26.cc:602:53: note: ‘dealii::DoFHandler<1>’ is not derived from ‘const dealii::hp::MappingCollection’ VectorTools::create_boundary_right_hand_side(dof_handler, ^In file included from /home/muhammad/software/dealii-8.5.0/examples/step-26/step-26.cc:46:0:/usr/local/include/deal.II/numerics/vector_tools.h:2136:8: note: candidate: template void dealii::VectorTools::create_boundary_right_hand_side(const dealii::hp::DoFHandler&, const dealii::hp::QCollection<(dim - 1)>&, const dealii::Function&, VectorType&, const std::set&) void create_boundary_right_hand_side ^/usr/local/include/deal.II/numerics/vector_tools.h:2136:8: note: template argument deduction/substitution failed:/home/muhammad/software/dealii-8.5.0/examples/step-26/step-26.cc:602:53: note: ‘dealii::DoFHandler<1>’ is not derived from ‘const dealii::hp::DoFHandler’ VectorTools::create_boundary_right_hand_side(dof_handler, ^make[3]: *** [CMakeFiles/step-26.dir/step-26.cc.o] Error 1CMakeFiles/step-26.dir/build.make:62: recipe for target 'CMakeFiles/step-26.dir/step-26.cc.o' failedmake[2]: *** [CMakeFiles/step-26.dir/all] Error 2make[3]: Leaving directory '/home/muhammad/software/dealii-8
Re: [deal.II] Inserting FESolution into systemmatrix
On 4/25/19 7:56 AM, Gabriel Peters wrote: > I wanted solution_o to be vector-valued. > And I want to multiply it with a gradient function. So I want > > local_matrix(i,j) += solution_o * fe_values[velocities].shape_gradient(i,q) > * > fe_values[velocities].shape_gradient(j,q); > > Does this look any better? Much :-) And if you take solution_o from the vector of tensors I mentioned, this should indeed work out. 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] Inserting FESolution into systemmatrix
Gabriel Peters Endenicher Str. 310 53121 Bonn 00491525/5478185 gabriel.pet...@koeln.de Am 25.04.19 um 15:28 schrieb Wolfgang Bangerth > On 4/25/19 4:51 AM, gabriel.pet...@koeln.de wrote: > > > > std::vector > solution_old_values(n_q_points); > > > > While assembling the system over all cells I call > > > > fe_values.get_function_values (solution_old, solution_old_values); > > > > In the Step-21 tutorial, they call > > > > const double old_s = old_solution_values[q](dim+1); > > > > to get the dim+2 component of the vector. > > > > Up to here, everything compiles fine, but I dont know how to apply the > > values > > correctly > > Therefore I have the following two questions: > > > > Do I simply have to call: const double solution_o = > > solution_old_values[q](0); > > to get the velocity component? In fact they are vector-valued, so a double > > doesn't make sense that much. > > Correct. The vector should be a `Tensor<1,dim>` and you would have to copy > the > first `dim` components of `solution_old_values[q]` into the tensor. That's > burdensome, so I would suggest that instead you use extractors: >std::vector> old_velocity_values(n_q_points); >fe_values[velocities].get_function_values (...); > and then you already have (only) the velocities and in their right data type. > > > > And thus my second question is: Do I apply the value correctly to the local > > matrix > > by calling: > > local_matrix(i,j) += solution_o * fe_values[velocities].shape_value(i,q) * > > fe_values[velocities].shape_value(j,q); > > This looks funny. You are multiplying a scalar 'solution_o' and phi_i*phi_j > where the latter two are vectors. This works, of course, but didn't you want > solution_o to be a vector as well? Oh thanks I Wrote bullshit. I wanted solution_o to be vector-valued. And I want to multiply it with a gradient function. So I want local_matrix(i,j) += solution_o * fe_values[velocities].shape_gradient(i,q) * fe_values[velocities].shape_gradient(j,q); Does this look any better? > > You are also missing the JxW. I left it out to here to keep Notation as Short as possible Best regards Gabriel > > 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. -- 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] Setting an initial condition by fixing the DOF values themselves
On 4/25/19 6:58 AM, Bruno Blais wrote: > I would like to be able to set-up an initial condition by fixing the value of > the DOF themselves by using their x,y(,z) position with a ParsedFuction. > > Generally, I set the initial condition by using an L2 projection of the > function (for complex function and higher order element I know this is more > appropriate). > However, for some stuff I would like to be able to fix the value of the DOF > directly. > My issue is, through the doxygen, I have not found a way to loop through the > DOF and to extract their position, even more when I have high order element > or > when I have a multiplicity (say velocity-vector + pressure). > > What I would be looking into doing would be something similar to : > > loop over local DOF; > get DOF position > get DOF component (u[0], u[1], u[2], p) > calculate parsed function using the point > set the value of the DOF > end > > Is it something that is possible? This is a weird way of iterating since I > know we always loop over the cells. Maybe there is a way to achieve this by > looping over the cells and then looping over the DOF of the cells? This would > be slightly redundant, but I do not care since I am doing this only once. You can get the locations of DoFs using DoFTools::map_dofs_to_support_points(). But maybe the easier way is to use VectorTools::interpolate_boundary_values() or project() and just pass it a function object that describes the function you want to set in terms of (x,y,z). This could be your ParsedFunction. interpolate_b_v() will then call your parsed function with a bunch of points that happen to correspond to the DoF locations, but it makes the connection internally and you no longer need to worry about which DoF sits where. 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] Inserting FESolution into systemmatrix
On 4/25/19 4:51 AM, gabriel.pet...@koeln.de wrote: > > std::vector > solution_old_values(n_q_points); > > While assembling the system over all cells I call > > fe_values.get_function_values (solution_old, solution_old_values); > > In the Step-21 tutorial, they call > > const double old_s = old_solution_values[q](dim+1); > > to get the dim+2 component of the vector. > > Up to here, everything compiles fine, but I dont know how to apply the values > correctly > Therefore I have the following two questions: > > Do I simply have to call: const double solution_o = solution_old_values[q](0); > to get the velocity component? In fact they are vector-valued, so a double > doesn't make sense that much. Correct. The vector should be a `Tensor<1,dim>` and you would have to copy the first `dim` components of `solution_old_values[q]` into the tensor. That's burdensome, so I would suggest that instead you use extractors: std::vector> old_velocity_values(n_q_points); fe_values[velocities].get_function_values (...); and then you already have (only) the velocities and in their right data type. > And thus my second question is: Do I apply the value correctly to the local > matrix > by calling: > local_matrix(i,j) += solution_o * fe_values[velocities].shape_value(i,q) * > fe_values[velocities].shape_value(j,q); This looks funny. You are multiplying a scalar 'solution_o' and phi_i*phi_j where the latter two are vectors. This works, of course, but didn't you want solution_o to be a vector as well? You are also missing the JxW. 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] nonlinear laplace problem with only dirichlet boundaries
On 4/25/19 3:20 AM, MD wrote: > I am trying to solve an nonlinear laplace problem with only dirichlet > boundaries using dealii with distributed triangulation. I do have a solution > vectorlocally_owned_sol and a newton_update. At the Dirichlet boundaries, I > want torestrict the newton update to be zero and the boundary values will be > taken in the locally_owned_sol. How can one achieved this in dealii. The > first > step is easy using the VectorTools::interpolate_boundary_values. I am > worried > about the second step how to construct it as well asinitalize the > newton_update with non_zero values at the Dirichlet boundaries. I don't think this is going to be any different than what is done in a sequential program. Have you taken a look at step-15? 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] Setting an initial condition by fixing the DOF values themselves
Hello everyone, I believe this is a relatively dumb question, but I seem to struggle with this basic concept. I would like to be able to set-up an initial condition by fixing the value of the DOF themselves by using their x,y(,z) position with a ParsedFuction. Generally, I set the initial condition by using an L2 projection of the function (for complex function and higher order element I know this is more appropriate). However, for some stuff I would like to be able to fix the value of the DOF directly. My issue is, through the doxygen, I have not found a way to loop through the DOF and to extract their position, even more when I have high order element or when I have a multiplicity (say velocity-vector + pressure). What I would be looking into doing would be something similar to : loop over local DOF; get DOF position get DOF component (u[0], u[1], u[2], p) calculate parsed function using the point set the value of the DOF end Is it something that is possible? This is a weird way of iterating since I know we always loop over the cells. Maybe there is a way to achieve this by looping over the cells and then looping over the DOF of the cells? This would be slightly redundant, but I do not care since I am doing this only once. Thank you for everything! Bruno -- 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: Error at Trilinos installation
It looks as if you have to recompile Trilinos using the -fPIC-switch to create position independent code. Am Donnerstag, 25. April 2019 12:10:52 UTC+2 schrieb gabrie...@koeln.de: > > Hey everyone, > I have an error at compiling deal_ii_with_Trilinos=ON, which I dont > understand. > > First I installed Trilinos,then I called > > cmake -DCMAKE_INSTALL_PREFIIX=/home/peters/Documents/installs > -DDEAL_II_WITH_TRILINOS=ON > -DTRILINOS_DIR=/home/peters/Documents/build_trilinos > > This all worked finely, but when I called make install I get the following > error: > > [ 55%] Linking CXX shared library ../lib/libdeal_II.so > /usr/bin/ld: > /home/peters/Documents/build_trilinos/lib/libmuelu-adapters.a(MueLu_CreateEpetraPreconditioner.cpp.o): > > relocation R_X86_64_PC32 against symbol > `_ZTVN7Teuchos17SerialTriDiMatrixIidEE' can not be used when making a > shared object; recompile with -fPIC > /usr/bin/ld: final link failed: Bad value > collect2: error: ld returned 1 exit status > source/CMakeFiles/deal_II.dir/build.make:964: recipe for target > 'lib/libdeal_II.so.9.1.0-pre' failed > make[2]: *** [lib/libdeal_II.so.9.1.0-pre] Error 1 > CMakeFiles/Makefile2:941: recipe for target > 'source/CMakeFiles/deal_II.dir/all' failed > make[1]: *** [source/CMakeFiles/deal_II.dir/all] Error 2 > Makefile:129: recipe for target 'all' failed > make: *** [all] Error 2 > peters@stud-04:~/programs/dealii$ m > > > Can somebody tell me, what went wrong? And where do I have to add the > order "recompile with -fPIC". > > Thanks a lot and best regards > > > Gabriel > -- 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] Inserting FESolution into systemmatrix
Hey everyone, I want to set up a BlockMatrix, with entries, which are enforced by a velocity components of aBlockVector "solution_old". The Vector contains velocity entries in the first dim terms and pressure entries in the last component. Similar to the 'assemble_system' function of Step-21 tutorial, I use std::vector > solution_old_values(n_q_points); While assembling the system over all cells I call fe_values.get_function_values (solution_old, solution_old_values); In the Step-21 tutorial, they call const double old_s = old_solution_values[q](dim+1); to get the dim+2 component of the vector. Up to here, everything compiles fine, but I dont know how to apply the values correctly Therefore I have the following two questions: Do I simply have to call: const double solution_o = solution_old_values[q](0); to get the velocity component? In fact they are vector-valued, so a double doesn't make sense that much. And thus my second question is: Do I apply the value correctly to the local matrix by calling: local_matrix(i,j) += solution_o * fe_values[velocities].shape_value(i,q) * fe_values[velocities].shape_value(j,q); As always, thanks a lot anf best regards Gabriel -- 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] Error at Trilinos installation
Hey everyone, I have an error at compiling deal_ii_with_Trilinos=ON, which I dont understand. First I installed Trilinos,then I called cmake -DCMAKE_INSTALL_PREFIIX=/home/peters/Documents/installs -DDEAL_II_WITH_TRILINOS=ON -DTRILINOS_DIR=/home/peters/Documents/build_trilinos This all worked finely, but when I called make install I get the following error: [ 55%] Linking CXX shared library ../lib/libdeal_II.so /usr/bin/ld: /home/peters/Documents/build_trilinos/lib/libmuelu-adapters.a(MueLu_CreateEpetraPreconditioner.cpp.o): relocation R_X86_64_PC32 against symbol `_ZTVN7Teuchos17SerialTriDiMatrixIidEE' can not be used when making a shared object; recompile with -fPIC /usr/bin/ld: final link failed: Bad value collect2: error: ld returned 1 exit status source/CMakeFiles/deal_II.dir/build.make:964: recipe for target 'lib/libdeal_II.so.9.1.0-pre' failed make[2]: *** [lib/libdeal_II.so.9.1.0-pre] Error 1 CMakeFiles/Makefile2:941: recipe for target 'source/CMakeFiles/deal_II.dir/all' failed make[1]: *** [source/CMakeFiles/deal_II.dir/all] Error 2 Makefile:129: recipe for target 'all' failed make: *** [all] Error 2 peters@stud-04:~/programs/dealii$ m Can somebody tell me, what went wrong? And where do I have to add the order "recompile with -fPIC". Thanks a lot and best regards Gabriel -- 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] Error at Trilinos installation
Hey everyone, I have an error at compiling deal_ii_with_Trilinos=ON, which I dont understand. First I installed,then I called cmake -DCMAKE_INSTALL_PREFIIX=/home/peters/Documents/installs -DDEAL_II_WITH_TRILINOS=ON -DTRILINOS_DIR=/home/peters/Documents/build_trilinos This all worked finely, but when I called make install I get the following error: [ 55%] Linking CXX shared library ../lib/libdeal_II.so /usr/bin/ld: /home/peters/Documents/build_trilinos/lib/libmuelu-adapters.a(MueLu_CreateEpetraPreconditioner.cpp.o): relocation R_X86_64_PC32 against symbol `_ZTVN7Teuchos17SerialTriDiMatrixIidEE' can not be used when making a shared object; recompile with -fPIC /usr/bin/ld: final link failed: Bad value collect2: error: ld returned 1 exit status source/CMakeFiles/deal_II.dir/build.make:964: recipe for target 'lib/libdeal_II.so.9.1.0-pre' failed make[2]: *** [lib/libdeal_II.so.9.1.0-pre] Error 1 CMakeFiles/Makefile2:941: recipe for target 'source/CMakeFiles/deal_II.dir/all' failed make[1]: *** [source/CMakeFiles/deal_II.dir/all] Error 2 Makefile:129: recipe for target 'all' failed make: *** [all] Error 2 peters@stud-04:~/programs/dealii$ m Can somebody tell me, what went wrong? And where do I have to add the order "recompile with -fPIC". Thanks a lot and best regards Gabriel -- 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] nonlinear laplace problem with only dirichlet boundaries
Hello all, I am trying to solve an nonlinear laplace problem with only dirichlet boundaries using dealii with distributed triangulation. I do have a solution vector locally_owned_sol and a newton_update. At the Dirichlet boundaries, I want to restrict the newton update to be zero and the boundary values will be taken in the locally_owned_sol. How can one achieved this in dealii. The first step is easy using the VectorTools::interpolate_boundary_values. I am worried about the second step how to construct it as well as initalize the newton_update with non_zero values at the Dirichlet boundaries. -- 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.
