[deal.II] initial velocity
Hello all,
This code just solves the wave equation. Despite the fact that I have
provided an initial velocity in line 106, it does not appear to be reading
it in lines 1024-1062. My solution looks like the wave equation with 0
initial velocity, i.e. one hump turns into two moving away from each other.
How can I ensure it is properly reading in the initial velocity I have
provided?
Thank you
--
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/e8535ad1-030c-4aea-b28d-bea7d83bf4dbn%40googlegroups.com.
/* -
*
* This file aims to solve a GR problem utilizing a Reissner-Nordstrom metric
*
* - */
#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 "helper.h"
#include
#include
using namespace dealii;
/* -- */
/* -- */
class Coefficients : public ParameterAcceptor
{
public:
typedef std::complex value_type;
static constexpr value_type imag{0., 1.};
Coefficients()
: ParameterAcceptor("A - Coefficients")
{
R_0 = 1.9;
add_parameter("R_0", R_0, "inner radius of the computational domain");
R_1 = 20.0;
add_parameter("R_1", R_1, "outer radius of the computational domain");
M = 1.0;
add_parameter("M", M, "mass of the black hole");
Q = 0.0;
add_parameter("Q", Q, "total charge of the black hole");
q = 1.;
add_parameter("q", q, "charge density");
A_param = 1.0;
add_parameter("A_param", A_param, "height of initial gaussian");
b_param = 10.0;
add_parameter("b_param", b_param, "position of peak of initial gaussian");
c_param = 2.0;
add_parameter("c_param", c_param, "wave speed of initial gaussian");
d_param = 1.0;
add_parameter("d_param", d_param, "standard deviation of initial gaussian");
f_param = 20.0;
add_parameter("f_param",f_param, "width of initial gaussian");
lam = 1.;
add_parameter("lambda", lam, "coefficient on nonlinear rhs.");
f = [&](double r) {
return 1. - (2. * M / r) + (Q * Q) / (r * r) + 0. * imag;
};
f_prime = [&](double r) {
return ((-2. * M)/ (r * r) + (2.*Q*Q)/(r*r*r)) + 0. * imag;
};
a = [&](double r) { return f(r)-2.; };
b = [&](double r) { return f_prime(r) + (2.*f(r)+2.*imag*q*Q-2.)/r; };
c = [&](double r) { return (2.*(1.-f(r))); };
d = [&](double r) { return (2.*f(r)+2.*imag*q*Q)/r; };
e = [&](double r) { return (imag*q*Q/(r*r)); };
/* IC = Aexp[-f{((x-b)+cT)^2/d^2}], initial_values = IC(0)*/
initial_values = [&](double r) {
//return -2.*(1.-20.*r + 2.*r*r)*std::exp(-(r-10)*(r-10))/(r*r*r); //For Manufactured Solution Run
return A_param*std::exp(-f_param * (r-b_param)*(r-b_param)/(pow(d_param,2.)));
};
initial_velocity = [&](double r) {
//return -4.*std::exp(-(r-10)*(r-10))/(r*r)-2.*(r-10)*initial_values(r);//For Manufactured Solution Run
return (-2.0 * c_param * f_param) / pow(d_param, 2.) * (r-b_param) * initial_values(r);
};
boundary_values = [&](double r, double t) {
//return 0.;
if (r == R_0)
return 0.;
if (r == R_1)
return 0.;
/* else
return 0.;*/
//return 0.;
//return .02 * std::exp(4. * ((1./(R_1-.9))-(1/(R_1+5.))+(4./(.9+5.;
//return 10. * std::exp(-(R_1-1.5)*(R_1-1.5)/2.);
};
MMS_rhs = [&](double r, double t){
return (2.+std::exp(-(r-10.+t)*(r-10.+t))*(1.+2.*(t-10.)*r+2.*r*r))/(r*r*r);
};
}
/* Publicly readable parameters: */
double R_0;
double R_1;
double lam;
double A_param;
double b_param;
double c_param;
double d_param;
double f_param;
/* Publicly readable functions: */
std::function a;
std::function b;
std::function c;
std::function d;
std::function e;
std::function f;
std::function initial_values;
std::function initial_velocity;
std::function boundary_values;
std::function MMS_rhs;
private:
/* Private functions: */
//std::function f;
std::function f_prime;
/* Parameters: */
double M;
double Q;
double q;
};
/* -- */
/* --
Re: [deal.II] boundaries in dealii
Hi, Thank you for taking a look. You are correct in that if I comment out VectorTools::interpolate_boundary_values in line 320, the solution does not change. If I understand you correctly, you are saying that the function apply_boundary_values defined in line 553 and called in line 630 must take in the system_matrix and must modify system_matrix in the same way I do the solution vector? This also means that when I call apply_boundary_values, I must move the location of this call until say line 646 or after I define system_matrix but before I solve? Thank you. On Tuesday, April 27, 2021 at 11:59:14 PM UTC-4 Wolfgang Bangerth wrote: > > Nikki, > apologies for not getting to this any earlier! > > > I've recently implemented a 1d standard wave equation and have been > messing > > around with different boundary conditions. For instance, I have run my > code > > with 1 dirichlet and 1 neumann, 2 neumann, and 1 absorbing condition and > 1 > > neumann (as in step-24). When I run step-23 from the examples, which > utilizes > > Dirichlet conditions, you can easily see the expected behavior of a > phase flip > > after interaction from the boundary, i.e. as the source approaches the > > boundary it looks like a cap and after the boundary interaction it looks > like > > a cup. However, in all 3 cases that I have tried, listed above, I am not > > observing the phase flip. No matter what boundary I change the left one > to > > (dirichlet, neumann, absorbing) it seems to be giving me the same > neumann > > behavior. Do you have any idea why this might be? > > I don't see any specific thing that looks wrong, but then it's an 800 line > code -- not easy to see anything without spending a couple of hours with > it. > > But here's where I would start: Look at the difference between Dirichlet > and > Neumann cases. If I see this right, setting up Dirichlet conditions > happens in > line 319 for boundary zero, which in line 162 you seem to think of as the > right boundary. So does the call to > VectorTools::interpolate_boundary_values > actually do anything? If you comment it out, does it change anything? > > I suspect no: You put the result into the 'boundary_value_map' variable, > but > that's a local variable whose lifetime ends at the end of the function, > and it > is never used after calling VectorTools::interpolate_boundary_values. So > that > looks like a bug -- or at least a pointless piece of code. > > Of course, you do the same again in line 573. There, you seem to set the > elements of the solution vector to the correct boundary values, but you do > not > touch the matrix. But you need to adjust the matrix before you solve with > it > if you want to apply the correct set of boundary values -- take a look at > step-3, step-4, step-26, for example. In particular the last one. > > > 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/518ad080-4b59-401c-9110-ebb4dd54dd57n%40googlegroups.com.
[deal.II] boundaries in dealii
I've recently implemented a 1d standard wave equation and have been messing
around with different boundary conditions. For instance, I have run my code
with 1 dirichlet and 1 neumann, 2 neumann, and 1 absorbing condition and 1
neumann (as in step-24). When I run step-23 from the examples, which
utilizes Dirichlet conditions, you can easily see the expected behavior of
a phase flip after interaction from the boundary, i.e. as the source
approaches the boundary it looks like a cap and after the boundary
interaction it looks like a cup. However, in all 3 cases that I have tried,
listed above, I am not observing the phase flip. No matter what boundary I
change the left one to (dirichlet, neumann, absorbing) it seems to be
giving me the same neumann behavior. Do you have any idea why this might be?
I have attached my code. If you were to run the code currently, it would
run with a dirichlet condition on the left and a neumann on the right. In
order to run with the absorbing boundary matrix , I comment out lines
316-320 and line 611 as well as uncommenting the very end of lines 616 and
645 which adds the boundary matrix to the computation.
Thank you!!
--
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/eb188eba-ba35-4324-ac20-7bee8db4d810n%40googlegroups.com.
/* -
*
* This file aims to solve the standard wave equation u_tt-u_rr=0 and variants
*
* - */
#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 "helper.h"
#include
#include
using namespace dealii;
/* -- */
/* -- */
class Coefficients : public ParameterAcceptor
{
public:
typedef std::complex value_type;
static constexpr value_type imag{0., 1.};
Coefficients()
: ParameterAcceptor("A - Coefficients")
{
R_0 = 0.0;
add_parameter("R_0", R_0, "inner radius of the computational domain");
R_1 = 10.0;
add_parameter("R_1", R_1, "outer radius of the computational domain");
A_param = 1.0;
add_parameter("A_param", A_param, "height of initial gaussian");
b_param = 5.0;
add_parameter("b_param", b_param, "position of peak of initial gaussian");
c_param = 2.0;
add_parameter("c_param", c_param, "wave speed of initial gaussian");
d_param = 1.0;
add_parameter("d_param", d_param, "standard deviation of initial gaussian");
f_param = 10.0;
add_parameter("f_param", f_param, "multiplicative constant of initial gaussian");
/* IC = A exp[-f{((x-b)+cT)^2/d^2}], initial_values1 = IC(0) */
initial_values1 = [&](double r) {
return A_param*std::exp(-f_param * (r-b_param)*(r-b_param)/(pow(d_param,2.)));
};
//IC_T = d/dT(IC)(0)
IC_T = [&](double r) {
return -(2.0 * c_param * f_param) / pow(d_param, 2.) * (r-b_param) * initial_values1(r);
};
boundary_values = [&](double r, double t) {
if (r == R_0)
return 0.;
//return t;
if (r == R_1)
return 0.;
};
}
/* Publicly readable parameters: */
double R_0;
double R_1;
double A_param;
double b_param;
double c_param;
double d_param;
double f_param;
/* Publicly readable functions: */
std::function initial_values1;
std::function IC_T;
std::function boundary_values;
};
/* -- */
/* -- */
template
class Discretization : public ParameterAcceptor
{
public:
static_assert(dim == 1, "Only implemented for 1D");
Discretization(const Coefficients &coefficients)
: ParameterAcceptor("B - Discretization")
, p_coefficients(&coefficients)
{
ParameterAcceptor::parse_parameters_call_back.connect(
std::bind(&Discretization::parse_parameters_callback, this));
refinement = 5;
add_parameter(
"refinement", refinement, "refinement of the spatial geometry");
order_mapping = 1;
add_parameter("order mapping", order_mapping, "order of the mapping");
order_finite_element = 1;
add_parameter("order finite element",
order_fi
Re: [deal.II] boundary matrix
I am using dealii version 9.1.1. The alternative you suggested has resolved
my issue.
Thank you,
Nikki
On Tuesday, March 9, 2021 at 4:16:36 AM UTC-5 Wolfgang Bangerth wrote:
> On 3/7/21 10:16 PM, Nikki Holtzer wrote:
> > Thank you for your input! That was my first attempt at tackling this.
> However,
> > I run into problems with the line
> >
> > for(const auto & face : cell-> face_iterators()).
> >
> > No matter what member function I choose off of the DoFCellAccessor class
> > template reference, I get some version of the following:
> >
> > class dealii::DoFCellAccessor, false>’ has no
> member
> > named ‘face_iterators’; did you mean ‘face_iterator’?
>
> I'm pretty sure that function should exist. What deal.II version are you
> using?
>
> But, if that doesn't work, you can always do this instead
> ```
> for(const auto & cell : triangulation.active_cell_iterators())
> for (unsigned int f=0; f::faces_per_cell; ++f)
> {
> const auto face = cell->face(f);
> if(face->at_boundary() && face->boundary_id() == 1)
> // do something useful
> }
> ```
>
> 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/7fa93ae5-ffab-4c25-96c7-3a81382c4fdan%40googlegroups.com.
Re: [deal.II] boundary matrix
Thank you for your input! That was my first attempt at tackling this.
However, I run into problems with the line
for(const auto & face : cell-> face_iterators()).
No matter what member function I choose off of the DoFCellAccessor class
template reference, I get some version of the following:
class dealii::DoFCellAccessor, false>’ has no
member named ‘face_iterators’; did you mean ‘face_iterator’?
If I then try the correction: I obtain the following:
invalid use of ‘using face_iterator'= class
dealii::TriaIterator,
false> >’
Like I said, this seems to happen with all the member functions I choose
for DoFCellAccessor.
On Sunday, March 7, 2021 at 3:07:09 PM UTC-5 peterr...@gmail.com wrote:
> I am not completely sure what you want to accomplish, but I think the
> following code should work for your purpose:
> ```
> for(const auto & cell : triangulation.active_cell_iterators())
> for(const auto & face : cell-> face_iterators())
> if(face->at_boundary() && face->boundary_id() == 1)
> // do something useful
> ```
> (plus minus a few typos)
>
> Hope that helps,
> Peter
> On Sunday, 7 March 2021 at 21:00:12 UTC+1 nhol...@math.arizona.edu wrote:
>
>> Yes, I have tried it. The problem I am having is with the
>> 'face->boundary_id() == 1' piece. I get an error which states that :
>> 'face' was not declared in this scope.
>>
>> The only other place that 'face' appears in my code is where I set the
>> boundary ids like the following:
>>
>> triangulation.begin_active->face()->set_boundary_id()
>>
>> Is there an easier way to access that boundary in the second half of that
>> if statement above other than
>>
>> (face->boundary_id() == 1)
>>
>> or is this the correct approach and it needs to be amended somehow?
>> On Sunday, March 7, 2021 at 1:49:12 PM UTC-5 Wolfgang Bangerth wrote:
>>
>>> On 3/7/21 5:30 AM, Nikki Holtzer wrote:
>>> >
>>> > For a 1d problem, I would like to implement a "boundary matrix" for an
>>> > absorbing boundary condition such as the one presented in Step 24 for
>>> only the
>>> > left boundary. Is it necessary to have the 'if' statement be a
>>> conditional
>>> > over faces like in Step 24 even though it is a 1d problem or is
>>> something like
>>> >
>>> > if (cell->at_boundary() && (face->boundary_id() == 1)) {
>>> > }
>>> > acceptable?
>>>
>>> The statement is ok I believe (have you tried?), the question is what
>>> you want
>>> to do after that.
>>>
>>> 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/69b308e3-13ef-43fb-a9d3-a765b7017fb0n%40googlegroups.com.
Re: [deal.II] boundary matrix
Yes, I have tried it. The problem I am having is with the
'face->boundary_id() == 1' piece. I get an error which states that :
'face' was not declared in this scope.
The only other place that 'face' appears in my code is where I set the
boundary ids like the following:
triangulation.begin_active->face()->set_boundary_id()
Is there an easier way to access that boundary in the second half of that
if statement above other than
(face->boundary_id() == 1)
or is this the correct approach and it needs to be amended somehow?
On Sunday, March 7, 2021 at 1:49:12 PM UTC-5 Wolfgang Bangerth wrote:
> On 3/7/21 5:30 AM, Nikki Holtzer wrote:
> >
> > For a 1d problem, I would like to implement a "boundary matrix" for an
> > absorbing boundary condition such as the one presented in Step 24 for
> only the
> > left boundary. Is it necessary to have the 'if' statement be a
> conditional
> > over faces like in Step 24 even though it is a 1d problem or is
> something like
> >
> > if (cell->at_boundary() && (face->boundary_id() == 1)) {
> > }
> > acceptable?
>
> The statement is ok I believe (have you tried?), the question is what you
> want
> to do after that.
>
> 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/f177b4d9-98a5-4c69-8a4e-e8e7ffe01827n%40googlegroups.com.
[deal.II] boundary matrix
Hello everyone!
For a 1d problem, I would like to implement a "boundary matrix" for an
absorbing boundary condition such as the one presented in Step 24 for only
the left boundary. Is it necessary to have the 'if' statement be a
conditional over faces like in Step 24 even though it is a 1d problem or is
something like
if (cell->at_boundary() && (face->boundary_id() == 1)) {
}
acceptable?
Thank you
--
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/c1cff077-e61f-4096-a75f-53c778deb38cn%40googlegroups.com.
Re: [deal.II] Linear Operator
Matthias, When you helped me implement the operator N_c which is added to my linear_operator , i.e. auto N_c = temp1_coef*M + temp_coef*Kronecker; auto system_matrix = linear_operator(linear_part) + N_c; the first term is not complete. Really I want, auto N_c = temp1_coef*|Psi^m|^2 *M + temp_coef*Kronecker; auto system_matrix = linear_operator(linear_part) + N_c; So I was trying in some way to represent that so I can add that term to 'temp_coef*Kronecker' which is of type linear operator. Nikki On Monday, October 26, 2020 at 3:07:39 PM UTC-4 Matthias Maier wrote: > Dear Nikki, > > A "linear operator" in deal.II is an object that linearly transforms a > vector into another one, i.e., it has a vector space as domain of > definition and a vector space as range space. > > In your case M.matrix_scalar_product(u,u) is (u^T M u), this is a > quadratic operation, not a linear one. You cannot represent this > operation with a linear operator. Did you intend to simply scale with > that value? > > Best, > Matthias > > > On Mon, Oct 26, 2020, at 13:40 CDT, Nikki Holtzer < > nhol...@math.arizona.edu> wrote: > > > Hello all, > > > > I am struggling to turn a matrix scalar product into a linear operator. > I > > currently have: > > > > temp = M.matrix_scalar_product(u,u); > > > > > > and I need to turn 'temp' into a linear operator. > > > > I have tried turning 'M' into a linear operator first and then > performing > > the multiplication and have also tried something like > > linear_operator(temp). Neither has compiled. > > > > Thank you in advance! > -- 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/f2d5c61a-49c2-4a7a-a984-1f3917b215bdn%40googlegroups.com.
[deal.II] Linear Operator
Hello all, I am struggling to turn a matrix scalar product into a linear operator. I currently have: temp = M.matrix_scalar_product(u,u); and I need to turn 'temp' into a linear operator. I have tried turning 'M' into a linear operator first and then performing the multiplication and have also tried something like linear_operator(temp). Neither has compiled. Thank you in advance! -- 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/2e306fd4-a862-4981-99ac-5e1aebb3ed63n%40googlegroups.com.
Re: [deal.II] FullMatrix to SparseMatrix
Yes it appears to be the system_matrix_inverse because I can print out the residual and it is indeed not empty. On Saturday, October 17, 2020 at 2:32:02 PM UTC-4 Wolfgang Bangerth wrote: > On 10/16/20 10:35 AM, Nikki Holtzer wrote: > > > > The error I'm receiving is: > > > > An error occurred in line <179> of file > > > > > > in function > > > > void dealii::FullMatrix::vmult(dealii::Vector&, > const > > dealii::Vector&, bool) const [with number2 = double; number > = double] > > > > The violated condition was: > > > > !this->empty() > > > > Additional information: > > > > (none) > > The error message states that you are trying to multiply a vector by a > matrix > that is empty (that is, has zero columns and/or rows). Have you run the > code > in a debugger to see which matrix this is? > > 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/5757b61c-20e4-4933-9376-d2301477a93cn%40googlegroups.com.
Re: [deal.II] FullMatrix to SparseMatrix
I was able to diagnose and fix the problem. Thank you. I still have one
remaining error since changing all my matrices to FullMatrix. I have
diagnosed that the issue is the very last line of code below, when
attempting to update my system. I have printed the line in red.
for (unsigned int m = 0; m < nonlinear_solver_limit; ++m) {
Vector residual = M_u * (new_solution - old_solution) + kappa *
(1. - theta) * S_u * new_solution
+ theta * kappa * S_u * old_solution - kappa * (1.- theta) *
lam * new_solution.norm_sqr() * new_solution - kappa * theta * lam *
old_solution.norm_sqr() * old_solution;
affine_constraints.set_zero(residual);
if (residual.linfty_norm() < nonlinear_solver_tol)
break;
double new_solution_normsq = -new_solution.norm_sqr();
const unsigned int k = new_solution.size();
FullMatrixKronecker_matrix(k,k);
for(unsigned int i=0; inonlinear_part_full;
nonlinear_part_full.copy_from(nonlinear_part);
FullMatrixnonlinear_part_inverse;
double nonlinear_part_coef = - kappa * (1.- theta) * lam *
new_solution_normsq;
double Kronecker_matrix_coef = -2.* lam * kappa * (1- theta);
nonlinear_part_full *= nonlinear_part_coef; // -k(1- theta) lambda
|Psi|^2*X
Kronecker_matrix *= Kronecker_matrix_coef; //(-2.* lam * kappa * (1-
theta)* Kronecker_matrix)
for(unsigned int i=0; i solver(solver_control);
auto system_matrix_inverse = inverse_operator(system_matrix, solver,
nonlinear_part_inverse);
Vector update = system_matrix_inverse * (-1. * residual);
The error I'm receiving is:
An error occurred in line <179> of file
in function
void dealii::FullMatrix::vmult(dealii::Vector&, const
dealii::Vector&, bool) const [with number2 = double; number =
double]
The violated condition was:
!this->empty()
Additional information:
(none)
I have also tried this approach
Vector update;
system_matrix_inverse.vmult(update, -1*residual)
but it resulted in the same error.
On Friday, October 16, 2020 at 10:37:40 AM UTC-4 Wolfgang Bangerth wrote:
>
> Nikki,
> I can't tell from these pieces of code what dimension is wrong. But you
> will
> find it very useful to run your program in a debugger so that you can see
> which matrix access it is that triggers the error, and why. There are a
> number
> of video lectures about running a program in a debugger.
>
> Best
> WB
>
>
> On 10/14/20 6:14 PM, Nikki Holtzer wrote:
> > The first error seems to be occurring at:
> >
> > affine_constraints.distribute_local_to_global(cell_nodal_mass_matrix,
> > local_dof_indices, nodal_mass_matrix);
> >
> > I have verified that cell_nodal_mass_matrix is indeed the same size as
> > nodal_mass_matrix.
> >
> > cell_nodal_mass_matrix is initialized as:
> >
> >
> > FullMatrix cell_nodal_mass_matrix(dofs_per_cell, dofs_per_cell);
> >
> >
> > and nodal_mass_matrix has the same dimensions.
> >
> > local_dof_indices seem to be the problem. I see that on dealii 's
> website that
> > as long as cell_nodal_mass_matrix and local_dof_indices have the same
> number
> > of elements this should work.
> >
> > I have defined local_dof_indices in the following way:
> >
> >
> > std::vector local_dof_indices(dofs_per_cell);
> >
> > fe_values.reinit(cell);
> >
> > cell->get_dof_indices(local_dof_indices);
> >
> >
> > The number of elements are:
> >
> > cell_nodal_mass_matrix: 136
> >
> > local_dof_indices: 24
> >
> > How would I resolve this?
> >
> > Thank you,
> >
> > Nikki
> >
> > On Wednesday, October 14, 2020 at 10:14:34 AM UTC-4 Nikki Holtzer wrote:
> >
> > Ok thank you that has resolved most of the adding issue. I am still
> > receiving the following error. It reads to me that there's still an error
> > on the 'FullMatrix::add()' function. I have verified that they are all
> > gone however. What else should I be looking for?
> >
> >
> >
> >
> > An error occurred in line <1335> of file
> > in
> function
> >
> > void dealii::FullMatrix<
> > >::add(dealii::FullMatrix< >::size_type,
> > dealii::FullMatrix< >::size_type, const
> > index_type*, const number2*, bool, bool) [with number2 = double;
> > index_type = unsigned int; number = double; dealii::FullMatrix<
> > >::size_type = long unsigned int]
> >
> > The violated condition was:
> >
> > static_cast > std::common_typem())>::type)>:: type>(row)
> > < sta
Re: [deal.II] FullMatrix to SparseMatrix
The first error seems to be occurring at: affine_constraints.distribute_local_to_global(cell_nodal_mass_matrix, local_dof_indices, nodal_mass_matrix); I have verified that cell_nodal_mass_matrix is indeed the same size as nodal_mass_matrix. cell_nodal_mass_matrix is initialized as: FullMatrix cell_nodal_mass_matrix(dofs_per_cell, dofs_per_cell); and nodal_mass_matrix has the same dimensions. local_dof_indices seem to be the problem. I see that on dealii 's website that as long as cell_nodal_mass_matrix and local_dof_indices have the same number of elements this should work. I have defined local_dof_indices in the following way: std::vector local_dof_indices(dofs_per_cell); fe_values.reinit(cell); cell->get_dof_indices(local_dof_indices); The number of elements are: cell_nodal_mass_matrix: 136 local_dof_indices: 24 How would I resolve this? Thank you, Nikki On Wednesday, October 14, 2020 at 10:14:34 AM UTC-4 Nikki Holtzer wrote: > Ok thank you that has resolved most of the adding issue. I am still > receiving the following error. It reads to me that there's still an error > on the 'FullMatrix::add()' function. I have verified that they are all gone > however. What else should I be looking for? > > > > > An error occurred in line <1335> of file > in > function > > void dealii::FullMatrix< > >::add(dealii::FullMatrix< >::size_type, > dealii::FullMatrix< >::size_type, const > index_type*, const number2*, bool, bool) [with number2 = double; index_type > = unsigned int; number = double; dealii::FullMatrix< > >::size_type = long unsigned int] > > The violated condition was: > > static_cast std::common_typem())>::type)>:: type>(row) < > static_cast std::common_typem())>::type)>:: > type>(this->m()) > > Additional information: > > Index 0 is not in the half-open range [0,0). In the current case, this > half-open range is in fact empty, suggesting that you are accessing an > element of an empty collection such as a vector that has not been set to > the correct size. > > > Stacktrace: > > --- > > #0 /home/u7/nholtzer/dealii-9.1.1/lib/libdeal_II.g.so.9.1.1: void > dealii::FullMatrix::add(unsigned long, > unsigned long, unsigned int const*, double const*, bool, bool) > > #1 /home/u7/nholtzer/dealii-9.1.1/lib/libdeal_II.g.so.9.1.1: void > dealii::AffineConstraints::distribute_local_to_global, > > dealii::Vector >(dealii::FullMatrix const&, > dealii::Vector const&, std::vector std::allocator > const&, dealii::FullMatrix&, > dealii::Vector&, bool, std::integral_constant) const > > #2 ./gravwave: void > dealii::AffineConstraints::distribute_local_to_global > > >(dealii::FullMatrix const&, std::vector std::allocator > const&, dealii::FullMatrix&) const > > #3 ./gravwave: OfflineData<1>::assemble() > > #4 ./gravwave: OfflineData<1>::prepare() > > #5 ./gravwave: main > > > On Tuesday, October 13, 2020 at 3:37:03 PM UTC-4 Wolfgang Bangerth wrote: > >> On 10/13/20 10:28 AM, Nikki Holtzer wrote: >> > The issue that I have when all matrices are FullMatrix matrices is I >> can't >> > seem to get the .add function to result in anything other than error. >> If I >> > keep ONLY both nonlinear_part and linear_part as FullMatrix >> > I receive: >> > >> > error: no matching function for call to >> > ‘dealii::FullMatrix::add(double, const >> dealii::SparseMatrix&)’ >> > >> >linear_part.add((1. - theta) * kappa, S_c); >> > >> > >> > Which leads me to believe that I cannot add a SparseMatrix (S_c) to a >> > FullMatrix (linear_part). >> >> Correct -- there is no such function. >> >> >> > I tried to back track and change S_c to be a >> > FullMatrix instead of a SparseMatrix >> > >> > For instance: >> > >> > >> > FullMatrix mass_matrix_unconstrainted; >> > mass_matrix_unconstrained = 0; >> > FullMatrix cell_mass_matrix(dofs_per_cell, dofs_per_cell); >> > cell_mass_matrix = 0; >> > >> > >> > Then I fill cell_mass_matrix. >> > >> > mass_matrix_unconstrained.add(local_dof_indices, cell_mass_matrix); >> > >> > error: no matching function for call to >> > ‘dealii::FullMatrix::add(std::vector&, >> > dealii::FullMatrix&)’ >> > >> > mass_matrix_unconstrained.add(local_dof_indices, >> cell_mass_matrix); >> >> Ri
Re: [deal.II] FullMatrix to SparseMatrix
Ok thank you that has resolved most of the adding issue. I am still receiving the following error. It reads to me that there's still an error on the 'FullMatrix::add()' function. I have verified that they are all gone however. What else should I be looking for? An error occurred in line <1335> of file in function void dealii::FullMatrix< >::add(dealii::FullMatrix< >::size_type, dealii::FullMatrix< >::size_type, const index_type*, const number2*, bool, bool) [with number2 = double; index_type = unsigned int; number = double; dealii::FullMatrix< >::size_type = long unsigned int] The violated condition was: static_castm())>::type)>:: type>(row) < static_castm())>::type)>:: type>(this->m()) Additional information: Index 0 is not in the half-open range [0,0). In the current case, this half-open range is in fact empty, suggesting that you are accessing an element of an empty collection such as a vector that has not been set to the correct size. Stacktrace: --- #0 /home/u7/nholtzer/dealii-9.1.1/lib/libdeal_II.g.so.9.1.1: void dealii::FullMatrix::add(unsigned long, unsigned long, unsigned int const*, double const*, bool, bool) #1 /home/u7/nholtzer/dealii-9.1.1/lib/libdeal_II.g.so.9.1.1: void dealii::AffineConstraints::distribute_local_to_global, dealii::Vector >(dealii::FullMatrix const&, dealii::Vector const&, std::vector > const&, dealii::FullMatrix&, dealii::Vector&, bool, std::integral_constant) const #2 ./gravwave: void dealii::AffineConstraints::distribute_local_to_global >(dealii::FullMatrix const&, std::vector > const&, dealii::FullMatrix&) const #3 ./gravwave: OfflineData<1>::assemble() #4 ./gravwave: OfflineData<1>::prepare() #5 ./gravwave: main On Tuesday, October 13, 2020 at 3:37:03 PM UTC-4 Wolfgang Bangerth wrote: > On 10/13/20 10:28 AM, Nikki Holtzer wrote: > > The issue that I have when all matrices are FullMatrix matrices is I > can't > > seem to get the .add function to result in anything other than error. If > I > > keep ONLY both nonlinear_part and linear_part as FullMatrix > > I receive: > > > > error: no matching function for call to > > ‘dealii::FullMatrix::add(double, const > dealii::SparseMatrix&)’ > > > >linear_part.add((1. - theta) * kappa, S_c); > > > > > > Which leads me to believe that I cannot add a SparseMatrix (S_c) to a > > FullMatrix (linear_part). > > Correct -- there is no such function. > > > > I tried to back track and change S_c to be a > > FullMatrix instead of a SparseMatrix > > > > For instance: > > > > > > FullMatrix mass_matrix_unconstrainted; > > mass_matrix_unconstrained = 0; > > FullMatrix cell_mass_matrix(dofs_per_cell, dofs_per_cell); > > cell_mass_matrix = 0; > > > > > > Then I fill cell_mass_matrix. > > > > mass_matrix_unconstrained.add(local_dof_indices, cell_mass_matrix); > > > > error: no matching function for call to > > ‘dealii::FullMatrix::add(std::vector&, > > dealii::FullMatrix&)’ > > > > mass_matrix_unconstrained.add(local_dof_indices, cell_mass_matrix); > > Right, that function also doesn't exist, but you could do > for (unsigned int i=0; i for (unsigned int j=0; j mass_matrix_unconstrainted(local_dof_indices[i], local_dof_indices[j]) > += cell_mass_matrix(i,j); > > > > Alternatively if I try: > > > > mass_matrix_unconstrained = mass_matrix_unconstrained + > > local_dof_indices*cell_mass_matrix; > > > > > > I receive. > > > > error: no match for ‘operator*’ (operand types are ‘std::vector int>’ > > and ‘dealii::FullMatrix’) > > Right. There is also no operator* for vector and > FullMatrix. It's also unclear to me what exactly that would > actually > supposed to do. > > Writing the double-loop above is the way to go. > > 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/54a09c9d-cfe6-4cd9-9899-f3138b1221c2n%40googlegroups.com.
Re: [deal.II] FullMatrix to SparseMatrix
The issue that I have when all matrices are FullMatrix matrices is I can't seem to get the .add function to result in anything other than error. If I keep ONLY both nonlinear_part and linear_part as FullMatrix I receive: error: no matching function for call to ‘dealii::FullMatrix::add(double, const dealii::SparseMatrix&)’ linear_part.add((1. - theta) * kappa, S_c); Which leads me to believe that I cannot add a SparseMatrix (S_c) to a FullMatrix (linear_part). I tried to back track and change S_c to be a FullMatrix instead of a SparseMatrix For instance: FullMatrix mass_matrix_unconstrainted; mass_matrix_unconstrained = 0; FullMatrix cell_mass_matrix(dofs_per_cell, dofs_per_cell); cell_mass_matrix = 0; Then I fill cell_mass_matrix. mass_matrix_unconstrained.add(local_dof_indices, cell_mass_matrix); error: no matching function for call to ‘dealii::FullMatrix::add(std::vector&, dealii::FullMatrix&)’ mass_matrix_unconstrained.add(local_dof_indices, cell_mass_matrix); Alternatively if I try: mass_matrix_unconstrained = mass_matrix_unconstrained + local_dof_indices*cell_mass_matrix; I receive. error: no match for ‘operator*’ (operand types are ‘std::vector’ and ‘dealii::FullMatrix’) mass_matrix_unconstrained = mass_matrix_unconstrained + local_dof_indices * cell_mass_matrix; I thought maybe this is a type error, i.e. I can't multiply a vector of int by a matrix of double but I can't seem to be able to convert the vector without error either. Thank you, Nikki On Monday, October 12, 2020 at 7:27:35 PM UTC-4 Wolfgang Bangerth wrote: > On 10/11/20 6:09 PM, Nikki Holtzer wrote: > > I was hoping to take the outer product which is formed in a FullMatrix, > called > > Kronecker_matrix below, multiply it by a coefficient, and add a matrix > to it > > called nonlinear_part_full. This entire operation I have tried to put > into a > > Sparse matrix, called Nonlinear_part, because I need to add this term to > > linear_part which is Sparse. I have done this simply because it was the > only > > way I could make the addition between linear_part (Sparse) and > > nonlinear_part_full (FullMatrix) work all while trying to make the code > as > > efficient as possible. > > But if the nonlinear part is a dense matrix, you will necessarily add > entries > into the sparse matrix that are not present in the sparsity pattern. When > you > build a sparse matrix with a certain sparsity pattern, you are saying > "these > are the entries that could possibly be nonzero, please allocate memory for > them". So when you later add the dense matrix to the sparse one, you're > breaking this promise. > > I don't quite understand why you copy everything into a sparse matrix in > your > code. Couldn't you just keep everything in dense matrices and do the > linear > solve with that? > > 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/a9ee9a8a-6cab-47b0-9ef1-49a133a6828dn%40googlegroups.com.
Re: [deal.II] FullMatrix to SparseMatrix
I have copied some code below.
I was hoping to take the outer product which is formed in a FullMatrix,
called Kronecker_matrix below, multiply it by a coefficient, and add a
matrix to it called nonlinear_part_full. This entire operation I have tried
to put into a Sparse matrix, called Nonlinear_part, because I need to add
this term to linear_part which is Sparse. I have done this simply because
it was the only way I could make the addition between linear_part (Sparse)
and nonlinear_part_full (FullMatrix) work all while trying to make the code
as efficient as possible.
template
void TimeStep::prepare()
{
std::cout << "TimeStep::prepare()" << std::endl;
const auto &offline_data = *p_offline_data;
linear_part.reinit(offline_data.sparsity_pattern);
auto &M_c = offline_data.mass_matrix;
auto &S_c = offline_data.stiffness_matrix;
linear_part.copy_from(M_c);
linear_part.add((1. - theta) * kappa, S_c);
linear_part_inverse.initialize(linear_part);
nonlinear_part.reinit(offline_data.sparsity_pattern);
auto &X = offline_data.inner_product_matrix;
nonlinear_part.copy_from(X);
}
template
void TimeStep::step(Vector &old_solution,double new_t)
{
const auto &offline_data = *p_offline_data;
const auto &coefficients = *offline_data.p_coefficients;
const auto &affine_constraints = offline_data.affine_constraints;
auto &manufactured = *p_manufactured;
const auto lam = coefficients.lam;
std::cout<<"lambda in Time step = "<> vector_memory;
typename VectorMemory>::Pointer
p_new_solution(vector_memory);
auto &new_solution = *p_new_solution;
new_solution = old_solution;
apply_boundary_values(offline_data, coefficients, new_t, new_solution,
old_solution);
for (unsigned int m = 0; m < nonlinear_solver_limit; ++m)
{ Vector residual = M_u * (new_solution - old_solution)
+ kappa * (1. - theta) * S_u * new_solution + theta * kappa *
S_u * old_solution - kappa * (1.- theta) * lam * new_solution.norm_sqr()
* new_solution - kappa * theta * lam *
old_solution.norm_sqr() * old_solution;
affine_constraints.set_zero(residual);
if (residual.linfty_norm() < nonlinear_solver_tol)
break;
double old_solution_normsq = -old_solution.norm_sqr();
const unsigned int k = old_solution.size();
FullMatrixKronecker_matrix(k,k);
for(unsigned int i=0; inonlinear_part_full;
nonlinear_part_full.copy_from(nonlinear_part);
double nonlinear_part_coef = - kappa * (1.- theta) * lam *
old_solution_normsq;
double Kronecker_matrix_coef = -2.* lam * kappa * (1- theta);
nonlinear_part_full *= nonlinear_part_coef; // -k(1- theta) lambda
|Psi|^2*X
Kronecker_matrix *= Kronecker_matrix_coef; //(-2.* lam * kappa * (1-
theta)* Kronecker_matrix)
nonlinear_part_full.add(1.,Kronecker_matrix); //-k(1- theta) lambda
|Psi|^2*X - 2k(1-theta)lambda Kronecker_matrix
SparseMatrix Nonlinear_part;
Nonlinear_part.reinit(offline_data.sparsity_pattern);
Nonlinear_part.copy_from(nonlinear_part_full);
linear_part.add(1., Nonlinear_part);
auto system_matrix = linear_operator(linear_part);
SolverControl solver_control(linear_solver_limit, linear_solver_tol);
SolverGMRES<> solver(solver_control);
auto system_matrix_inverse = inverse_operator(system_matrix, solver,
linear_part_inverse);
Vector update = system_matrix_inverse * (-1. * residual);
affine_constraints.set_zero(update);
new_solution += update;
}
On Sunday, October 11, 2020 at 5:38:12 PM UTC-4 Wolfgang Bangerth wrote:
>
> Nikki,
> we can continue to speculate, but your description of what is happening is
> not
> specific enough to really know. You'll have to show us the code you use to
> be
> sure.
>
> I don't understand the connection between the code you show and the outer
> product. The outer product is a dense matrix, but the sparsity pattern you
> create is sparse. Were you hoping to put the outer product of two vectors
> into
> a sparse matrix?
>
> Best
> W.
>
>
> On 10/10/20 3:27 PM, Nikki Holtzer wrote:
> > My sparsity pattern was created in the following way:
> >
> > DynamicSparsityPattern c_sparsity(dof_handler.n_dofs(),
> dof_handler.n_dofs());
> >
> > DoFTools::make_sparsity_pattern(dof_handler, c_sparsity,
> > affine_constraints, true);
> >
> >
> > which I have used throughout my code with no problems.
> >
> >
> > The FullMatrix, Mat1, that I form comes from the outer product of vec1
> and
> > vec2 both of size dof_handler.n_dofs.
> >
> > I create a Sparse Matrix, Mat2 that i initialize with the sparsity
> pattern
> > written above. I am not seeing
Re: [deal.II] FullMatrix to SparseMatrix
My sparsity pattern was created in the following way: DynamicSparsityPattern c_sparsity(dof_handler.n_dofs(), dof_handler.n_dofs()); DoFTools::make_sparsity_pattern(dof_handler, c_sparsity, affine_constraints, true); which I have used throughout my code with no problems. The FullMatrix, Mat1, that I form comes from the outer product of vec1 and vec2 both of size dof_handler.n_dofs. I create a Sparse Matrix, Mat2 that i initialize with the sparsity pattern written above. I am not seeing why Mat1 of size dof_handler.n_dofs X dof_handler.n_dofs does not fit into a sparsity pattern that's defined as DynamicSparsityPattern c_sparsity(dof_handler.n_dofs(), dof_handler.n_dofs()); I am indeed getting an error that there are not enough entries in the SparseMatrix, Mat2. Thank you! On Wednesday, October 7, 2020 at 4:58:45 PM UTC-4 Wolfgang Bangerth wrote: > On 10/7/20 2:51 PM, Nikki Holtzer wrote: > > > > When doing so I receive the following error: > > > > > > An error occurred in line <195> of file > > in function > > > > dealii::SparseMatrix& > > dealii::SparseMatrix::operator=(double) [with number = double] > > > > The violated condition was: > > > > cols != nullptr > > > > Additional information: > > > > (none) > > You haven't associated a sparsity pattern with the sparse matrix yet. > You'll > have to set that first, and you need to make sure that it has an entry for > each nonzero in the full matrix you are hoping to copy over! > > 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/56cca082-6ae2-4b51-af0d-85ed95661ac4n%40googlegroups.com.
Re: [deal.II] FullMatrix to SparseMatrix
Hello everyone, I was just curious if anyone had any suggestions for the above 'const' qualifier problem I have stated above? Thanks! Nikki On Thursday, October 8, 2020 at 9:16:24 AM UTC-4 Nikki Holtzer wrote: > When I do the above, with the modifications you have provided, I still > obtain the following: > > error: no match for ‘operator+’ (operand types are ‘const > dealii::SparseMatrix’ and ‘dealii::SparseMatrix’) > > or > > error: no matching function for call to > ‘dealii::SparseMatrix::add(double, dealii::SparseMatrix&) > const’ > > > I have not added the 'const' qualifier to my definitions of either of my > sparse matrices. > > On Wednesday, October 7, 2020 at 4:58:45 PM UTC-4 Wolfgang Bangerth wrote: > >> On 10/7/20 2:51 PM, Nikki Holtzer wrote: >> > >> > When doing so I receive the following error: >> > >> > >> > An error occurred in line <195> of file >> > in >> function >> > >> > dealii::SparseMatrix& >> > dealii::SparseMatrix::operator=(double) [with number = double] >> > >> > The violated condition was: >> > >> > cols != nullptr >> > >> > Additional information: >> > >> > (none) >> >> You haven't associated a sparsity pattern with the sparse matrix yet. >> You'll >> have to set that first, and you need to make sure that it has an entry >> for >> each nonzero in the full matrix you are hoping to copy over! >> >> 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/ad114115-9411-410f-829d-5d7fac91c24cn%40googlegroups.com.
Re: [deal.II] FullMatrix to SparseMatrix
When I do the above, with the modifications you have provided, I still obtain the following: error: no match for ‘operator+’ (operand types are ‘const dealii::SparseMatrix’ and ‘dealii::SparseMatrix’) or error: no matching function for call to ‘dealii::SparseMatrix::add(double, dealii::SparseMatrix&) const’ I have not added the 'const' qualifier to my definitions of either of my sparse matrices. On Wednesday, October 7, 2020 at 4:58:45 PM UTC-4 Wolfgang Bangerth wrote: > On 10/7/20 2:51 PM, Nikki Holtzer wrote: > > > > When doing so I receive the following error: > > > > > > An error occurred in line <195> of file > > in function > > > > dealii::SparseMatrix& > > dealii::SparseMatrix::operator=(double) [with number = double] > > > > The violated condition was: > > > > cols != nullptr > > > > Additional information: > > > > (none) > > You haven't associated a sparsity pattern with the sparse matrix yet. > You'll > have to set that first, and you need to make sure that it has an entry for > each nonzero in the full matrix you are hoping to copy over! > > 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/c597e3c0-a53e-444e-ba60-8db28efb2d4fn%40googlegroups.com.
[deal.II] FullMatrix to SparseMatrix
Hello everyone, I am trying to recast a current FullMatrix as a SparseMatrix in order to add this matrix to another SparseMatrix. I have tried something like this: SparseMatrix Mat_Sparse; Mat_Sparse.copy_from(Mat_Full); Mat_Sparse_2.add(1.,Mat_Sparse); When doing so I receive the following error: An error occurred in line <195> of file in function dealii::SparseMatrix& dealii::SparseMatrix::operator=(double) [with number = double] The violated condition was: cols != nullptr Additional information: (none) How should I amend this recasting? Thank you! -- 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/63583cc9-55bd-4b13-8689-3afda516192an%40googlegroups.com.
Re: [deal.II] outer product of two vectors
Thank you! On Tuesday, October 6, 2020 at 8:31:47 PM UTC-4 Wolfgang Bangerth wrote: > On 10/6/20 5:50 PM, Nikki Holtzer wrote: > > > > I am trying to form a cross product/ outer product of two vectors of > type > > deallii:Vector. I have attempted to use some of the built in > functions > > for the outer product from the Tensor Class but have had no luck. I > can't seem > > to get anything other than > > > > error: no matching function for call to 'outer_product(vec1, vec2);' > > > > I have tried recasting my vec1/vec2 as Tensors but have run into a > similar > > error message. > > > > Is there a built in vector cross product? Alternatively, how could I > recast my > > vectors and then use the built in functions from the Tensor Class and > finally > > recast them back into vectors? > > The easy way is to do > > const unsigned int n = vec.size(); > FullMatrix o_p (n,n); > for (unsigned int i=0; i for (unsigned int j=0; j o_p(i,j) = vec[i] * vec[j]; > > But the issue is that generally you end up with a full matrix this way. Is > that what you want? How large are your vectors? > > 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/a7420864-4151-43ba-b3f5-0ef23a12c371n%40googlegroups.com.
[deal.II] outer product of two vectors
Hello everyone, I am trying to form a cross product/ outer product of two vectors of type deallii:Vector. I have attempted to use some of the built in functions for the outer product from the Tensor Class but have had no luck. I can't seem to get anything other than error: no matching function for call to 'outer_product(vec1, vec2);' I have tried recasting my vec1/vec2 as Tensors but have run into a similar error message. Is there a built in vector cross product? Alternatively, how could I recast my vectors and then use the built in functions from the Tensor Class and finally recast them back into vectors? Thank you! -- 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/bc43013d-1159-4ebe-bd9b-7364e8afbba8n%40googlegroups.com.
