Re: [deal.II] Applying non-homogeneous Neumann BC in step-18
On 7/15/23 11:33, Mohammad Amir Kiani Fordoei wrote: const Tensor<1,dim> n = fe_face_values.normal_vector(q_point); cell_rhs (i)+=( (n*traction_component*fe_face_values.shape_value(i, q_point))* fe_face_values.JxW(q_point)); }//for q_point }//for dof_per_cell But in computation of cell rhs i get this error: /home/amir/eclipse-workspace/mech1/mech.cc:983:22: error: no match for ‘operator+=’ (operand types are ‘double’ and ‘dealii::Tensor<1, 3>’) 983 | cell_rhs (i)+=( (n*traction_component*fe_face_values.shape_value(i, q_point))* | ^~ 984 | fe_face_values.JxW(q_point)); As error says, the normal vector with 3 component cannot be multiplied by double values in other terms. No, that's not actually what the error message says :-) It doesn't say anything about multiplication. If you read it carefully, it says that you cannot assign a right hand side of type Tensor<1,3> to a left hand side of type double. That makes sense. The types you have are because your right hand side is the product of n -- a tensor traction_component -- not sure, but probably a scalar fe_face_values.shape_value(i, q_point)) -- a scalar, see the function's documentation fe_face_values.JxW(q_point) -- a scalar What you *want* is a vector-valued shape function. You get that by "extracting" this from the fe_values object in the same way as is discussed in the vector-valued documentation module. I think something like this should work: FEValuesExtractors::Vector displacement(0); cell_rhs (i)+=( (n*traction_component*fe_face_values[displacement].shape_value(i, q_point))* fe_face_values.JxW(q_point)); 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. To view this discussion on the web visit https://groups.google.com/d/msgid/dealii/9f48ee6b-cf61-499b-985e-268524743ee0%40colostate.edu.
Re: [deal.II] Integrate difference isn't picking on the
On 7/16/23 07:22, Abbas Ballout wrote: I get this error when I call the exact_solution_vector_function.gradient(p1): !ERROR Message: /An error occurred in line <135> of file in function dealii::Tensor<1, dim, Number> dealii::FunctionRangeNumberType>::gradient(const dealii::Point&, unsigned int) const [with int dim = 2; RangeNumberType = double] The violated condition was: false Additional information: You (or a place in the library) are trying to call a function that is declared as a virtual function in a base class but that has not been overridden in your derived class./ / / Calling value works fine but not gradient even though I have overridden both functions. You do in your own class, but VectorFunctionFromTensorFunction does not. That's not by design, but probably simply because nobody has ever had a need for this. Would you like to try and write a patch that adds VectorFunctionFromTensorFunction::gradient()? You could model it on VectorFunctionFromTensorFunction::value(). Best Wolfgang -- 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. To view this discussion on the web visit https://groups.google.com/d/msgid/dealii/528f591f-ba68-0a94-d8d7-948665e1e119%40colostate.edu.
Re: [deal.II] Integrate difference isn't picking on the
In hindsight I could've written something simpler. Here it is: #include #include using namespace dealii; template class Exact : public TensorFunction<1, dim> { public: Exact() : TensorFunction<1, dim>(dim) {} virtual Tensor<1, dim> value(const Point &p) const override; virtual Tensor<2, dim> gradient(const Point &p) const override; }; template Tensor<1, dim> Exact::value(const Point &p) const { return Tensor<1, dim>(); } template Tensor<2, dim> Exact::gradient(const Point &p) const { return Tensor<2, dim>(); } int main() { const int dim = 2; Exact exact_solution; VectorFunctionFromTensorFunction exact_solution_vector_function(exact_solution, 0, dim); Point p1; exact_solution_vector_function.value(p1); exact_solution_vector_function.gradient(p1); return 0; } I get this error when I call the exact_solution_vector_function.gradient (p1): !ERROR Message: *An error occurred in line <135> of file in functiondealii::Tensor<1, dim, Number> dealii::Function::gradient(const dealii::Point&, unsigned int) const [with int dim = 2; RangeNumberType = double]The violated condition was: falseAdditional information: You (or a place in the library) are trying to call a function that isdeclared as a virtual function in a base class but that has not beenoverridden in your derived class.* Calling value works fine but not gradient even though I have overridden both functions. On Friday, July 14, 2023 at 8:55:34 PM UTC+2 Wolfgang Bangerth wrote: > On 7/14/23 09:34, Abbas Ballout wrote: > > // THis throws exception > > VectorTools::integrate_difference(dof_handler, > > solution, > > exact_solution_vector_function, > > difference_per_cell, > > quadrature_formula, > > VectorTools::H1_seminorm);; > > Well, what's the error message? A good rule of thumb is that 50% debugging > consists of carefully reading what the error message says. > > 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. To view this discussion on the web visit https://groups.google.com/d/msgid/dealii/c1157723-d4a3-4c67-96f7-8eb98aa57fd5n%40googlegroups.com.