Re: [deal.II] Discrepancy in Computing the Inf-Sup Constants Using Matlab and deal.II

2024-06-03 Thread Najwa Alshehri
Wolfgang:

Thank you for bringing this to my attention. I would definitely love to 
share the code in the code gallery for the benefit of others. In this case, 
I will create a code for a simpler problem, such as the Stokes problem or 
mixed Laplace problem, to illustrate the concept of the inf-sup condition 
and how to compute the constant. As soon as it is ready, I will submit it 
there.

Best regards,
Najwa

On Sunday, June 2, 2024 at 5:54:23 PM UTC+3 Wolfgang Bangerth wrote:

>
> Najwa:
> on a separate note: Computing the inf-sup constant is an interesting 
> problem. 
> It would be nice to have a program that shows how to do that available as 
> the 
> base for others' experiments. Could I interest you in submitting your 
> code, 
> once you've got it working, to the code gallery?
> https://dealii.org/code-gallery.html
>
> Best
> W.
>

-- 
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/9db91592-1dcc-494d-9fa4-eadb0b3d130an%40googlegroups.com.


Re: [deal.II] Discrepancy in Computing the Inf-Sup Constants Using Matlab and deal.II

2024-06-02 Thread Najwa Alshehri
Dear Luca,

Thank you for your answer, that works perfectly. 

Best,
Najwa

On Sunday, June 2, 2024 at 12:27:54 PM UTC+3 luca@gmail.com wrote:

> To be more specific, you can use inverse_operator, to which you pass A, a 
> solver, and a preconditioner, and it returns an operator the applies the 
> inverse of A to a vector. Just as your mass inverse.
>
> Luca
>
> Il giorno 2 giu 2024, alle ore 11:25, Luca Heltai  ha 
> scritto:
>
> Dear Najwa, 
>
>
> You don’t need to create the inverse of A. Just its action, pretty much in 
> the exact same way you did for the mass matrix, but replacing the inverse 
> operator done with umfpack with one done with CG. 
>
> Luca
>
> Il giorno 2 giu 2024, alle ore 09:42, Najwa Alshehri  
> ha scritto:
>
> Wolfgang and all,
>
> I have a positive definite matrix M, and its inverse can be quickly and 
> easily found in deal.II using the following lines:
>
> M_inv_umfpack.initialize(M);
> auto M = linear_operator(M);
> auto M_inv = linear_operator(M, M_inv_umfpack);
>
> However, when solving the direct system, I need to provide the inverse of 
> AA to the Arpack solver. This requires solving the system using the 
> Conjugate Gradient (CG) method a number of times equal to the size of AA. 
> As the mesh is refined, this size increases, and depending on the finite 
> element method (FE) I am using, it could be even larger.
>
> While I have been able to obtain at least the first three results for many 
> FEs I am testing which is enough for the current study, I am seeking a 
> faster solution for a general eigenvalue problem of the same type. If 
> anyone has suggestions for a more efficient approach to finding the inverse 
> of AA, I would greatly appreciate it.
>
> One more time, you have always been great supporters. This group is very 
> useful. Thank you to you Wolfgang, to Mathias for the answer, and to all 
> members and developers of deal.ii.
>
> Best regards,
> Najwa
>
>
>
> On Saturday, June 1, 2024 at 11:57:45 PM UTC+3 Wolfgang Bangerth wrote:
>
>> On 6/1/24 14:49, Najwa Alshehri wrote: 
>> > 
>> > I decided to solve the exact problem directly, namely AA x = \lambda M 
>> x. To 
>> > achieve this, I computed the inverse of the matrix AA= Bt * A^inv *  B 
>> using a 
>> > Conjugate Gradient (CG) solver. Subsequently, I solved for the exact 
>> > eigenvalues, and to my satisfaction, I obtained results that aligned 
>> with 
>> > those obtained from MATLAB. 
>> > 
>> > If anyone has insights into why solving for the reciprocal yields 
>> different 
>> > results, I would greatly appreciate your input on this matter. 
>>
>> Najwa, 
>> I'm glad you figured it out because I don't really have a good idea about 
>> what 
>> ARPACK does. Typically, 
>> AA x = \lambda M x 
>> is solved in exactly this way because M may be a matrix that is singular. 
>> In 
>> many applications (not sure whether yours is one of those), M is a matrix 
>> that 
>> can have a large null-space -- for example, if you are computing 
>> eigenvalues 
>> of the Stokes operator, in some applications you have that M = [M_u 0 ; 0 
>> 0], 
>> that is, it has a large null block on all pressures. This matrix is not 
>> invertible, whereas AA is. By convention, the invertible matrix is kept 
>> on the 
>> left side of the eigenvalue equation. 
>>
>> 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+un...@googlegroups.com.
> To view this discussion on the web visit 
> https://groups.google.com/d/msgid/dealii/e448210c-be8d-4460-8b41-27a4425f7accn%40googlegroups.com
>  
> <https://groups.google.com/d/msgid/dealii/e448210c-be8d-4460-8b41-27a4425f7accn%40googlegroups.com?utm_medium=email_source=footer>
> .
>
>

-- 
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/fb1efd0c-f161-44f6-a863-99c560ae4886n%40googlegroups.com.


Re: [deal.II] Discrepancy in Computing the Inf-Sup Constants Using Matlab and deal.II

2024-06-02 Thread Najwa Alshehri
Wolfgang and all,

I have a positive definite matrix M, and its inverse can be quickly and 
easily found in deal.II using the following lines:

M_inv_umfpack.initialize(M);
auto M = linear_operator(M);
auto M_inv = linear_operator(M, M_inv_umfpack);

However, when solving the direct system, I need to provide the inverse of 
AA to the Arpack solver. This requires solving the system using the 
Conjugate Gradient (CG) method a number of times equal to the size of AA. 
As the mesh is refined, this size increases, and depending on the finite 
element method (FE) I am using, it could be even larger.

While I have been able to obtain at least the first three results for many 
FEs I am testing which is enough for the current study, I am seeking a 
faster solution for a general eigenvalue problem of the same type. If 
anyone has suggestions for a more efficient approach to finding the inverse 
of AA, I would greatly appreciate it.

One more time, you have always been great supporters. This group is very 
useful. Thank you to you Wolfgang, to Mathias for the answer, and to all 
members and developers of deal.ii.

Best regards,
Najwa



On Saturday, June 1, 2024 at 11:57:45 PM UTC+3 Wolfgang Bangerth wrote:

> On 6/1/24 14:49, Najwa Alshehri wrote:
> > 
> > I decided to solve the exact problem directly, namely AA x = \lambda M 
> x. To 
> > achieve this, I computed the inverse of the matrix AA= Bt * A^inv *  B 
> using a 
> > Conjugate Gradient (CG) solver. Subsequently, I solved for the exact 
> > eigenvalues, and to my satisfaction, I obtained results that aligned 
> with 
> > those obtained from MATLAB.
> > 
> > If anyone has insights into why solving for the reciprocal yields 
> different 
> > results, I would greatly appreciate your input on this matter.
>
> Najwa,
> I'm glad you figured it out because I don't really have a good idea about 
> what 
> ARPACK does. Typically,
> AA x = \lambda M x
> is solved in exactly this way because M may be a matrix that is singular. 
> In 
> many applications (not sure whether yours is one of those), M is a matrix 
> that 
> can have a large null-space -- for example, if you are computing 
> eigenvalues 
> of the Stokes operator, in some applications you have that M = [M_u 0 ; 0 
> 0],
> that is, it has a large null block on all pressures. This matrix is not 
> invertible, whereas AA is. By convention, the invertible matrix is kept on 
> the 
> left side of the eigenvalue equation.
>
> 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/e448210c-be8d-4460-8b41-27a4425f7accn%40googlegroups.com.


Re: [deal.II] Discrepancy in Computing the Inf-Sup Constants Using Matlab and deal.II

2024-06-01 Thread Najwa Alshehri
Hello all,

One more time, thanks for your support.  I've *managed to identify and 
resolve the issue I was facing*, and I'd like to share the solution here 
for the benefit of others who might encounter similar challenges.

The problem stemmed from my approach of solving for the reciprocals of the 
eigenvalues, which I did just because the inverse of the mass matrix M was 
easier. Specifically,  I was solving M x = 1/\lambda AA x. The question of why 
this was the problem remains unsolved.

I decided to solve the exact problem directly, namely AA x = \lambda M x.  To 
achieve this, I computed the inverse of the matrix AA= Bt * A^inv *  B  using 
a Conjugate Gradient (CG) solver. Subsequently, I solved for the exact 
eigenvalues, and to my satisfaction, I obtained results that aligned with 
those obtained from MATLAB.

If anyone has insights into why solving for the reciprocal yields different 
results, I would greatly appreciate your input on this matter.

Best,
Najwa
On Friday, May 31, 2024 at 11:41:27 AM UTC+3 Najwa Alshehri wrote:

> Dear Wolfgang:
>
> Thank you for your answer. I have checked the matrices and I printed the 
> case for a 2 by 2 mesh from both Matlab and deal.II. With a simple change 
> of rows and columns, they are exactly the same ( see the files attached). 
> Not sure yet why they result in different results of the eigenvalue 
> problem. There is a very importent detail that I want to point out here. 
> When I solved on the very course mesh (2 by 2) I got the following error, 
> which is not the case when I solved on a 4 by 4 mesh for example.
>
> An error occurred in line <791> of file 
>  in function
> void dealii::ArpackSolver::solve(const MatrixType1&, const 
> MatrixType2&, const INVERSE&, std::vector >&, 
> std::vector&, unsigned int) [with VectorType = 
> dealii::Vector; MatrixType1 = 
> dealii::LinearOperator, dealii::Vector, 
> dealii::internal::LinearOperatorImplementation::EmptyPayload>; MatrixType2 
> = dealii::LinearOperator, dealii::Vector, 
> dealii::internal::LinearOperatorImplementation::EmptyPayload>; INVERSE = 
> dealii::SparseDirectUMFPACK]
> The violated condition was: 
> false
> Additional information: 
> Error with dnaupd, info -2. Check documentation of ARPACK
>
>
> Dear Mathias:
>
> Thank you for your answer. The matrices are exactly the same in both 
> codes. The different is that I solve for the recipricals in deal.II. To be 
> precise,
>
> let AA= Bt A^inv B 
> and let M be a mass matrix
>
> where B, A, M are the same matrices in both Matlab and deal.II.
>
> In Matlab: I solve eig(AA, M) (AA x = \lambda M x)  then I seeks the 
>  square root of the smalles \lambda 
> In dealii: I solve using ArpackSolver solving M x = 1/\Lambda AA x) 
> so I do first 
> M_inv_umfpack.initialize(M);
> auto M = linear_operator(M);
> auto M_inv = linear_operator(M, M_inv_umfpack);
> then, I create a linear operator of the other matrices and the inverse of 
> A and perform AA= Bt * A^inv *  B 
> lastly I solve : arpack_solver.solve( M, AA, M_inv_umfpack, 
> eigenvalues,eigenvectors);
>
> So I am seeking the largest eigenvalue and \lambda = 1/largest eigenvalue 
> and the infsup constant is the square root of that.
>
> Could solving for the recipricals is the issue? if yes what should I fill 
> in the third intery of the solving if Iam doing 
> arpack_solver.solve( AA, M, , eigenvalues,eigenvectors);
>
> One more time, your questions and help are highly appreciated.
> Best,
> Najwa
> On Friday, May 31, 2024 at 12:40:04 AM UTC+3 Matthias Maier wrote:
>
>> Do you use the same "mass matrix" in both cases, meaning do you call 
>> eig(A, M) where A is the stiffness matrix and M is the mass matrix in 
>> matlab? 
>>
>>
>> You will need to compute the spectrum of the eigenvalue problem 
>>
>> A x = \lambda M x 
>>
>>
>> Best, 
>> Matthias 
>>
>>
>> On Thu, May 30, 2024, at 10:55 CDT, Najwa Alshehri  
>> wrote: 
>>
>> > Dear all, 
>> > 
>> > 
>> > Thank you for your always support and help. 
>> > 
>> > 
>> > I am encountering an issue while solving an eigenvalue problem related 
>> to 
>> > the computation of the discrete inf-sup constant for a saddle point 
>> > problem. Specifically, I am solving the following system: 
>> > 
>> > 
>> > Bt A^-1 B eigenvector = eigenvalue M eigenvector, 
>> > 
>> > 
>> > where Bt A^-1 B is clearly symmetric and at least semi-positive 
>> definite, 
>> > and M is a symmetric positive definite matrix. 
>> > 
>> > 
>> > The pr

Re: [deal.II] Discrepancy in Computing the Inf-Sup Constants Using Matlab and deal.II

2024-05-31 Thread Najwa Alshehri
Dear Wolfgang:

Thank you for your answer. I have checked the matrices and I printed the 
case for a 2 by 2 mesh from both Matlab and deal.II. With a simple change 
of rows and columns, they are exactly the same ( see the files attached). 
Not sure yet why they result in different results of the eigenvalue 
problem. There is a very importent detail that I want to point out here. 
When I solved on the very course mesh (2 by 2) I got the following error, 
which is not the case when I solved on a 4 by 4 mesh for example.

An error occurred in line <791> of file 
 in function
void dealii::ArpackSolver::solve(const MatrixType1&, const 
MatrixType2&, const INVERSE&, std::vector >&, 
std::vector&, unsigned int) [with VectorType = 
dealii::Vector; MatrixType1 = 
dealii::LinearOperator, dealii::Vector, 
dealii::internal::LinearOperatorImplementation::EmptyPayload>; MatrixType2 
= dealii::LinearOperator, dealii::Vector, 
dealii::internal::LinearOperatorImplementation::EmptyPayload>; INVERSE = 
dealii::SparseDirectUMFPACK]
The violated condition was: 
false
Additional information: 
Error with dnaupd, info -2. Check documentation of ARPACK


Dear Mathias:

Thank you for your answer. The matrices are exactly the same in both codes. 
The different is that I solve for the recipricals in deal.II. To be precise,

let AA= Bt A^inv B 
and let M be a mass matrix

where B, A, M are the same matrices in both Matlab and deal.II.

In Matlab: I solve eig(AA, M) (AA x = \lambda M x)  then I seeks the 
 square root of the smalles \lambda 
In dealii: I solve using ArpackSolver solving M x = 1/\Lambda AA x) 
so I do first 
M_inv_umfpack.initialize(M);
auto M = linear_operator(M);
auto M_inv = linear_operator(M, M_inv_umfpack);
then, I create a linear operator of the other matrices and the inverse of A 
and perform AA= Bt * A^inv *  B 
lastly I solve : arpack_solver.solve( M, AA, M_inv_umfpack, 
eigenvalues,eigenvectors);

So I am seeking the largest eigenvalue and \lambda = 1/largest eigenvalue 
and the infsup constant is the square root of that.

Could solving for the recipricals is the issue? if yes what should I fill 
in the third intery of the solving if Iam doing 
arpack_solver.solve( AA, M, , eigenvalues,eigenvectors);

One more time, your questions and help are highly appreciated.
Best,
Najwa
On Friday, May 31, 2024 at 12:40:04 AM UTC+3 Matthias Maier wrote:

> Do you use the same "mass matrix" in both cases, meaning do you call
> eig(A, M) where A is the stiffness matrix and M is the mass matrix in
> matlab?
>
>
> You will need to compute the spectrum of the eigenvalue problem
>
> A x = \lambda M x
>
>
> Best,
> Matthias
>
>
> On Thu, May 30, 2024, at 10:55 CDT, Najwa Alshehri  
> wrote:
>
> > Dear all,
> >
> >
> > Thank you for your always support and help.
> >
> >
> > I am encountering an issue while solving an eigenvalue problem related 
> to 
> > the computation of the discrete inf-sup constant for a saddle point 
> > problem. Specifically, I am solving the following system:
> >
> >
> > Bt A^-1 B eigenvector = eigenvalue M eigenvector,
> >
> >
> > where Bt A^-1 B is clearly symmetric and at least semi-positive 
> definite, 
> > and M is a symmetric positive definite matrix.
> >
> >
> > The problem arises when comparing the inf-sup constants obtained using 
> > Matlab and dealii (solving on the same mesh ):
> >
> >
> > 1. Matlab: Using the "eig" solver, the inf-sup constant decays to zero 
> > as the mesh is refined.
> > 2. deal.II: Using the "ArpackSolver", the inf-sup constant is bounded 
> > away from zero.
> >
> > The discrepancy is puzzling, as both methods should give consistent 
> > results. To illustrate this issue, I have attached a figure comparing 
> the 
> > two results and included in this post the two data files obtained from 
> > Matlab and deal.II.
> >
> >
> > I would greatly appreciate any insights or suggestions on what might be 
> > causing the discrepancy.
> >
> >
> > Your help is highly appreciated,
> >
> > Najwa
>

-- 
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/ed61d756-3f6c-4d1c-a7ae-7337def7ebc6n%40googlegroups.com.
Mu2 =

0.7229   -0.1385 0   -0.1385   -0.3193 0 0 
0

[deal.II] Discrepancy in Computing the Inf-Sup Constants Using Matlab and deal.II

2024-05-30 Thread Najwa Alshehri


Dear all,


Thank you for your always support and help.


I am encountering an issue while solving an eigenvalue problem related to 
the computation of the discrete inf-sup constant for a saddle point 
problem. Specifically, I am solving the following system:


Bt A^-1 B eigenvector = eigenvalue M eigenvector,


where Bt A^-1 B is clearly symmetric and at least semi-positive definite, 
and M is a symmetric positive definite matrix.


The problem arises when comparing the inf-sup constants obtained using 
Matlab and dealii (solving on the same mesh ):


   1. Matlab: Using the "eig" solver, the inf-sup constant decays to zero 
   as the mesh is refined.
   2. deal.II: Using the "ArpackSolver", the inf-sup constant is bounded 
   away from zero.

The discrepancy is puzzling, as both methods should give consistent 
results. To illustrate this issue, I have attached a figure comparing the 
two results and included in this post the two data files obtained from 
Matlab and deal.II.


I would greatly appreciate any insights or suggestions on what might be 
causing the discrepancy.


Your help is highly appreciated,

Najwa


-- 
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/d28832e9-9f39-4300-b8c6-ef0504f2bc38n%40googlegroups.com.
0.25 0.346321
0.125 0.215878
0.0625 0.167721
0.03125 0.153921
0.2500 0.0601
0.1250 0.0162
0.0625 0.0041
0.0312  0.0010


inf_sup.pdf
Description: Adobe PDF document


[deal.II] Issue with intersection of two meshes with curved boundaries

2023-11-24 Thread Najwa Alshehri


Dear developers and users,

 

I have two meshes one is immersed in the other. I wanted to find the 
intersection between the two meshes, so I used the following function.

 

NonMatching::compute_intersection(omega_grid_tools_cache,

  omega2_grid_tools_cache,

  4,  // degree

  1e-20); // to


include header 


source code 



This function uses CGAL to find a quadrature formula to integrate exactly 
on the intersection of two meshes which neglecte intersections of areas 
with tolerance smaller than “tol“ that one chooses and gives a quadrature 
formula on the triangulation of the intersection area that integrates 
exactly polynomials of a specific degree ( which allows maximum degree of 
4).


Say that I want to integrate, on the intersection area, a polynomial of 
order 4.

 

I noticed that If I am considering a circular immersed domain (unlike 
square or L-shaped domains), after a few cycles, the quadrature formula is 
not accurate enough. To be precise, I find that the sum of the weights of 
the quadrature formula defined on the triangulation of the exact 
intersection of the two meshes does not sum up (up to a tolerance) to the 
measure of the domain. When this occurs, the solution that is solved by 
evaluating the integral considering the quadrature formula on the exact 
intersection is no longer correct and the error starts to diverge in the 
later cycles after this point.

 

Moreover, the difference gets large suddenly, in one cycle, the difference 
was relatively smaller (1e-13), and in the next, it is much larger (1e-8) 
as can be seen in the attached plot (Plot shows the difference between the 
sum of the weights on the whole domain defined on the triangulation of the 
exact intersection of the two meshes and the measure of the domain under 
uniform refinement of the mesh). 

 

Any suggestion on what could be the issue and what should I do to fix it?

 

Thanks

Najwa

-- 
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/a5044577-482d-4ba8-9adb-7a5a379668a7n%40googlegroups.com.


test_quad.pdf
Description: Adobe PDF document


[deal.II] Re: Q: deserialize after updating dealii

2023-11-07 Thread Najwa Alshehri
Dear Bruno,

I appreciate your observation. I reprocessed the solution by serializing it 
with the updated version of Dealii, and subsequently, when I deserialized 
it in a different code, it functioned as expected. 

I appreciate your help.
Best,
Najwa

On Tuesday, November 7, 2023 at 10:49:18 AM UTC+3 Najwa Alshehri wrote:

> Hello Bruno,
>
> I have the solution serialized with the old version. Do you think this is 
> the issue? I will do everything with the new version today and update you 
> if it works.
>
> Best,
> Najwa
>
> On Monday, November 6, 2023 at 5:33:11 PM UTC+3 bruno.t...@gmail.com 
> wrote:
>
>> Najwa,
>>
>> Are you serializing/deserializing the solution using the same version of 
>> deal.II and Boost? If you serialized the solution using the old version of 
>> deal.II, I don't think you can deserialize it using a newer version of 
>> deal.II
>>
>> Best,
>>
>> Bruno
>>
>> On Monday, November 6, 2023 at 3:03:07 AM UTC-5 najwaa...@gmail.com 
>> wrote:
>>
>>> Dear developers and group members,
>>>
>>>
>>> I was using an old version of Dealii where I was serializing some 
>>> solutions and deserializing them where needed. In the code where I 
>>> deserialize the solution I used to use the following headers:
>>>
>>>
>>> #include 
>>> #include 
>>> #include 
>>>
>>> #include 
>>> #include 
>>>
>>> and I used 
>>> std::ifstream iss("serialized_ex2");
>>> boost::archive::text_iarchive ia(iss);
>>>
>>> It used to work fine. However, after updating Dealii, the code breaks 
>>> and I get the following error:
>>>
>>> Exception on processing: 
>>> unregistered class
>>> Aborting!
>>>
>>> Do you have an idea of what could be the missing header or class?
>>>
>>> Your help is highly appreciated,
>>> Najwa
>>>
>>>
>>>
>>>

-- 
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/c1713f83-e6e9-4e01-9fd4-f8d6f2330ca7n%40googlegroups.com.


[deal.II] Re: Q: deserialize after updating dealii

2023-11-06 Thread Najwa Alshehri
Hello Bruno,

I have the solution serialized with the old version. Do you think this is 
the issue? I will do everything with the new version today and update you 
if it works.

Best,
Najwa

On Monday, November 6, 2023 at 5:33:11 PM UTC+3 bruno.t...@gmail.com wrote:

> Najwa,
>
> Are you serializing/deserializing the solution using the same version of 
> deal.II and Boost? If you serialized the solution using the old version of 
> deal.II, I don't think you can deserialize it using a newer version of 
> deal.II
>
> Best,
>
> Bruno
>
> On Monday, November 6, 2023 at 3:03:07 AM UTC-5 najwaa...@gmail.com wrote:
>
>> Dear developers and group members,
>>
>>
>> I was using an old version of Dealii where I was serializing some 
>> solutions and deserializing them where needed. In the code where I 
>> deserialize the solution I used to use the following headers:
>>
>>
>> #include 
>> #include 
>> #include 
>>
>> #include 
>> #include 
>>
>> and I used 
>> std::ifstream iss("serialized_ex2");
>> boost::archive::text_iarchive ia(iss);
>>
>> It used to work fine. However, after updating Dealii, the code breaks and 
>> I get the following error:
>>
>> Exception on processing: 
>> unregistered class
>> Aborting!
>>
>> Do you have an idea of what could be the missing header or class?
>>
>> Your help is highly appreciated,
>> Najwa
>>
>>
>>
>>

-- 
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/059d4f62-a3bc-4120-859d-85f3bb1a8ac6n%40googlegroups.com.


[deal.II] Q: deserialize after updating dealii

2023-11-06 Thread Najwa Alshehri
Dear developers and group members,


I was using an old version of Dealii where I was serializing some solutions 
and deserializing them where needed. In the code where I deserialize the 
solution I used to use the following headers:


#include 
#include 
#include 

#include 
#include 

and I used 
std::ifstream iss("serialized_ex2");
boost::archive::text_iarchive ia(iss);

It used to work fine. However, after updating Dealii, the code breaks and I 
get the following error:

Exception on processing: 
unregistered class
Aborting!

Do you have an idea of what could be the missing header or class?

Your help is highly appreciated,
Najwa



-- 
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/83727530-49d8-421c-a0db-80ad4ab3519en%40googlegroups.com.


Re: [deal.II] Q: Question about extracting part of a vector

2023-08-03 Thread Najwa Alshehri


Dear team,

 

I want to express my gratitude for the support you have provided thus far. 
I would like to provide you with an update on the current situation and 
seek your assistance regarding a couple of issues.

 

Firstly, I have successfully implemented the missing function 
`VectorFunctionFromScalarFunctionObject::gradient()` in my project. 
Although it appears to be functioning correctly, the execution speed is 
extremely slow. In fact, it took until the next working day to obtain the 
results for only two cycles of refinement. Based on the response from 
@Wolfgang, it seems that the issue lies with the FeFieldFunction.

 

Allow me to provide some context. I have a reference solution serialized on 
a fine mesh. In my current code, my objective is to deserialize this 
solution and utilize it as the exact solution. To achieve this, I have 
employed the FeFieldFunction. Now, I have two scenarios in which I need to 
utilize the FeFieldFunction.

 

The first scenario involves directly using the FeFieldFunction when my 
computed solution consists of a single component. In this case, I have 
successfully calculated the L2 and H1 norms of the error using the 
`integrate_difference` and `compute_global_error` functions. The 
calculations proceeded smoothly, and the execution speed was reasonable.

 

However, the second scenario involves using the FeFieldFunction when the 
computed solution consists of two components. My objective is to calculate 
the error for the first component only. To accomplish this, I employed the 
new VectorFunctionFromScalarFunctionObject class locally. Unfortunately, 
the computation of the L2 norm in this case was considerably slow and the 
H1 norm was very slow.

 

Given the current situation, I have two questions that I hope you can 
assist me with, as you have done in the past:

 

1. Regarding my specific project, is there an alternative method, other 
than the FeFieldFunction, to save the deserialized solution as a function 
and utilize it as the exact solution? It is crucial that the solution be 
saved in the first component of a vector with two components.

 

2. In a more general sense, how can I test this function to determine if 
the sluggishness is solely attributable to the FeFieldFunction? 
Additionally, I would like to verify whether the function is functioning 
correctly or not, to be able to send it in a pull request to the dealii git 
repository.

 

Your guidance and expertise in resolving these issues would be greatly 
appreciated.

 

Thank you once again for your continued support.

 

Best regards,

Najwa

On Tuesday, August 1, 2023 at 12:50:27 PM UTC+3 abbas.b...@gmail.com wrote:

> Najwa, 
> The files I edited are in the pull request here: 
> https://github.com/dealii/dealii/pull/15805.
> Let me know if I can help in a better way or of there is something I can 
> do. 
>
> Abbas. 
>
> On Tuesday, August 1, 2023 at 5:48:41 AM UTC+2 najwaa...@gmail.com wrote:
>
>> Hello Abbas,
>>
>> Thank you for your response. I have actually written something, but I 
>> need to test it. It would be great to have a look at your work as well, so 
>> we can compare.
>>
>> Best regards,
>> Najwa
>>
>> On Mon, 31 Jul 2023 at 10:48 PM Abbas Ballout  
>> wrote:
>>
>>> Najwaa, 
>>>
>>> I submitted a pull request <https://github.com/dealii/dealii/pull/15805> 
>>> recently about something similar 
>>> <https://groups.google.com/g/dealii/c/H0kLG8RJKVc>. (I guess)
>>> Maybe that would help. 
>>>
>>> Abbas 
>>>  
>>>
>>> On Monday, July 31, 2023 at 7:37:59 PM UTC+2 najwaa...@gmail.com wrote:
>>>
>>>> Dear Wolfgang,
>>>>
>>>> Thank you for your answer. I will try to write a patch to the deal.II 
>>>> sources that implement the missing function. This might require some time. 
>>>> I will come back here if needed.
>>>>
>>>>
>>>> Thank you all for your quick answers,
>>>> Najwa
>>>>
>>>> On Wednesday, July 26, 2023 at 11:41:15 PM UTC+3 Wolfgang Bangerth 
>>>> wrote:
>>>>
>>>>> On 7/26/23 12:17, Najwa Alshehri wrote: 
>>>>> > 
>>>>> > Regarding the function "VectorFucntionFromScalarFunctionObject," I 
>>>>> have 
>>>>> > observed that it works for computing the L2 norm of the error. 
>>>>> However, when I 
>>>>> > attempted to compute the H1_seminorm of the error, I noticed that 
>>>>> the gradient 
>>>>> > was not implemented for this function. Please let me know if my 
>>>>> understanding 
>>>>> > is incorrect, as I referre

Re: [deal.II] Q: Question about extracting part of a vector

2023-07-31 Thread Najwa Alshehri
Hello Abbas,

Thank you for your response. I have actually written something, but I need
to test it. It would be great to have a look at your work as well, so we
can compare.

Best regards,
Najwa

On Mon, 31 Jul 2023 at 10:48 PM Abbas Ballout 
wrote:

> Najwaa,
>
> I submitted a pull request <https://github.com/dealii/dealii/pull/15805>
> recently about something similar
> <https://groups.google.com/g/dealii/c/H0kLG8RJKVc>. (I guess)
> Maybe that would help.
>
> Abbas
>
>
> On Monday, July 31, 2023 at 7:37:59 PM UTC+2 najwaa...@gmail.com wrote:
>
>> Dear Wolfgang,
>>
>> Thank you for your answer. I will try to write a patch to the deal.II
>> sources that implement the missing function. This might require some time.
>> I will come back here if needed.
>>
>>
>> Thank you all for your quick answers,
>> Najwa
>>
>> On Wednesday, July 26, 2023 at 11:41:15 PM UTC+3 Wolfgang Bangerth wrote:
>>
>>> On 7/26/23 12:17, Najwa Alshehri wrote:
>>> >
>>> > Regarding the function "VectorFucntionFromScalarFunctionObject," I
>>> have
>>> > observed that it works for computing the L2 norm of the error.
>>> However, when I
>>> > attempted to compute the H1_seminorm of the error, I noticed that the
>>> gradient
>>> > was not implemented for this function. Please let me know if my
>>> understanding
>>> > is incorrect, as I referred to the  for this.
>>>
>>> Yes, correct. Someone needs to implement
>>> VectorFucntionFromScalarFunctionObject::gradient()
>>> and/or
>>> VectorFucntionFromScalarFunctionObject::gradient_list()
>>> in exactly the same way as the value() and value_list() functions are
>>> implemented.
>>>
>>>
>>> > Could you kindly advise me on the simplest way to accomplish this?
>>>
>>> The easiest way is to create a class derived from
>>> VectorFucntionFromScalarFunctionObject in your own project that
>>> implements the
>>> missing function. The better way is to write a patch to the deal.II
>>> sources
>>> that implements the missing function because then others can use your
>>> work in
>>> the future as well :-)
>>>
>>>
>>> > On a side note, although "VectorFucntionFromScalarFunctionObject"
>>> works with
>>> > the L2 norm of the error, it is quite computationally intensive.
>>>
>>> It isn't VectorFucntionFromScalarFunctionObject that's expensive, but
>>> the
>>> FEFieldFunction. That's sort of inherently the case with this class. If
>>> you
>>> can, it should be avoided and you should try an evaluate the FE field
>>> you have
>>> at regular quadrature points, via FEValues.
>>>
>>> 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 a topic in the
> Google Groups "deal.II User Group" group.
> To unsubscribe from this topic, visit
> https://groups.google.com/d/topic/dealii/DqH2auneWaY/unsubscribe.
> To unsubscribe from this group and all its topics, send an email to
> dealii+unsubscr...@googlegroups.com.
> To view this discussion on the web visit
> https://groups.google.com/d/msgid/dealii/dd5f93f8-fb58-45c6-9ab2-91ef6c5bdf56n%40googlegroups.com
> <https://groups.google.com/d/msgid/dealii/dd5f93f8-fb58-45c6-9ab2-91ef6c5bdf56n%40googlegroups.com?utm_medium=email_source=footer>
> .
>

-- 
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/CABQM8dT6n4f2r-Efj4XR%3DegyzXPW9UtY3Z9_vxiX1MeUhKNLWg%40mail.gmail.com.


Re: [deal.II] Q: Question about extracting part of a vector

2023-07-31 Thread Najwa Alshehri
Dear Wolfgang,

Thank you for your answer. I will try to write a patch to the deal.II 
sources that implement the missing function. This might require some time. 
I will come back here if needed.


Thank you all for your quick answers,
Najwa

On Wednesday, July 26, 2023 at 11:41:15 PM UTC+3 Wolfgang Bangerth wrote:

> On 7/26/23 12:17, Najwa Alshehri wrote:
> > 
> > Regarding the function "VectorFucntionFromScalarFunctionObject," I have 
> > observed that it works for computing the L2 norm of the error. However, 
> when I 
> > attempted to compute the H1_seminorm of the error, I noticed that the 
> gradient 
> > was not implemented for this function. Please let me know if my 
> understanding 
> > is incorrect, as I referred to the  for this.
>
> Yes, correct. Someone needs to implement
> VectorFucntionFromScalarFunctionObject::gradient()
> and/or
> VectorFucntionFromScalarFunctionObject::gradient_list()
> in exactly the same way as the value() and value_list() functions are 
> implemented.
>
>
> > Could you kindly advise me on the simplest way to accomplish this?
>
> The easiest way is to create a class derived from 
> VectorFucntionFromScalarFunctionObject in your own project that implements 
> the 
> missing function. The better way is to write a patch to the deal.II 
> sources 
> that implements the missing function because then others can use your work 
> in 
> the future as well :-)
>
>
> > On a side note, although "VectorFucntionFromScalarFunctionObject" works 
> with 
> > the L2 norm of the error, it is quite computationally intensive.
>
> It isn't VectorFucntionFromScalarFunctionObject that's expensive, but the 
> FEFieldFunction. That's sort of inherently the case with this class. If 
> you 
> can, it should be avoided and you should try an evaluate the FE field you 
> have 
> at regular quadrature points, via FEValues.
>
> 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/e1aa5c5a-46cc-417b-89d5-c1bbfd47a641n%40googlegroups.com.


Re: [deal.II] Q: Question about extracting part of a vector

2023-07-26 Thread Najwa Alshehri
Greetings everyone,

I want to express my gratitude for your assistance with my question, which 
has evolved into something more profound.

Regarding the function "VectorFucntionFromScalarFunctionObject," I have 
observed that it works for computing the L2 norm of the error. However, 
when I attempted to compute the H1_seminorm of the error, I noticed that 
the gradient was not implemented for this function. Please let me know if 
my understanding is incorrect, as I referred to the 
 for this.

I am curious if there is a way to resolve this issue. My idea is to obtain 
the gradient of the "FeFeildFunction" since I know that the gradient is 
implemented for this function, as per 
. Afterward, I can define a new 
function within a class, based on "VectorFucntionFromScalarFunctionObject," 
where the gradient of that function would be the gradient of the 
"FeFeildFunction."

Could you kindly advise me on the simplest way to accomplish this? 

On a side note, although "VectorFucntionFromScalarFunctionObject" works 
with the L2 norm of the error, it is quite computationally intensive.

Best regards,
Najwa
On Sunday, July 23, 2023 at 4:47:50 PM UTC+3 Najwa Alshehri wrote:

> Dear Luca,
>
> This is interesting. I am exited to try this out. 
>
> Thank you so much,
> Najwa
>
> On Sunday, July 23, 2023 at 3:39:05 PM UTC+3 luca@gmail.com wrote:
>
>> FEFieldFunction has the same number of components of the underlying 
>> dofhandler. If that has one component, then so does the FEFieldFunction. 
>> What you could do is wrap the FEFieldFunction into another function, e.g.
>>
>> The deal.II Library: VectorFunctionFromScalarFunctionObject< dim, 
>> RangeNumberType > Class Template Reference 
>> <https://www.dealii.org/current/doxygen/deal.II/classVectorFunctionFromScalarFunctionObject.html>
>> dealii.org 
>> <https://www.dealii.org/current/doxygen/deal.II/classVectorFunctionFromScalarFunctionObject.html>
>> [image: deal.ico] 
>> <https://www.dealii.org/current/doxygen/deal.II/classVectorFunctionFromScalarFunctionObject.html>
>>  
>> <https://www.dealii.org/current/doxygen/deal.II/classVectorFunctionFromScalarFunctionObject.html>
>>
>> You could pass a lambda function as an argument to the constructor of the 
>> function, that uses your FEFieldFunction object, I.e.,
>>
>> VectorFucntionFromScalarFuncObj my_fun([&](const auto ){ return 
>> fe_field_func.value(p);}, 0, 2);
>>
>> Luca
>>
>> Il giorno 21 lug 2023, alle ore 22:05, Najwa Alshehri <
>> najwaa...@gmail.com> ha scritto:
>>
>> Dear Daniel, your answer makes sense, Finally it worked.
>>
>>
>> Dear Luca, thank you for your answer. Indeed, the ExactSolution is a 
>> class I built myself. I followed your suggestion and it was the magic line 
>> that solved the issue. *However,* say that I have in another code 
>> instead of using a ExactSolution that I build myself, I used an 
>> FeFieldFunction for a reference solution, how would I make this function 
>> allocate two components?
>>
>> Best,
>> Najwa
>>
>> On Friday, July 21, 2023 at 1:54:42 PM UTC+3 luca@gmail.com wrote:
>>
>>> Dear Najwa, 
>>>
>>> how many components do the `ExactSolution2` function has? 
>>>
>>> The error complains about `exact_solution.n_components, n_components` 
>>> not being equal, in particular it tells you that `ExactSolution2` only 
>>> has one component. 
>>>
>>> If you built the ExactSolution class yourself, make sure at construction 
>>> time you tell Function to allocate two components: 
>>>
>>> ExactSolution2::ExactSolution2( … ) 
>>> : Functions::Function(2) 
>>>
>>> Best, 
>>> Luca. 
>>>
>>>
>>> > On 20 Jul 2023, at 12:15, Najwa Alshehri  wrote: 
>>> > 
>>> > Hello again, 
>>> > 
>>> > I have a follow-up question. Does this ComponentSelectFunction work 
>>> also with vectors that are not written as blocked vectors? I have applied 
>>> it before when I was dealing with blocked vectors and it worked perfectly. 
>>> > 
>>> > So, I did this, 
>>> > const ComponentSelectFunction primal_mask(0,2); 
>>> > Later, 
>>> > VectorTools::integrate_difference(omega2_dh, 
>>> > u_omega2, 
>>> > ExactSolution2(), 
>>> > cellwise_errors2, 
>>> > quadrature, 
>>> > VectorTools::L2_norm 
>>> > ,_mask); 
>>> > const double u2_l2_error = 
>>> > Vec

Re: [deal.II] Q: Question about extracting part of a vector

2023-07-23 Thread Najwa Alshehri
Dear Luca,

This is interesting. I am exited to try this out. 

Thank you so much,
Najwa

On Sunday, July 23, 2023 at 3:39:05 PM UTC+3 luca@gmail.com wrote:

> FEFieldFunction has the same number of components of the underlying 
> dofhandler. If that has one component, then so does the FEFieldFunction. 
> What you could do is wrap the FEFieldFunction into another function, e.g.
>
> The deal.II Library: VectorFunctionFromScalarFunctionObject< dim, 
> RangeNumberType > Class Template Reference 
> <https://www.dealii.org/current/doxygen/deal.II/classVectorFunctionFromScalarFunctionObject.html>
> dealii.org 
> <https://www.dealii.org/current/doxygen/deal.II/classVectorFunctionFromScalarFunctionObject.html>
> [image: deal.ico] 
> <https://www.dealii.org/current/doxygen/deal.II/classVectorFunctionFromScalarFunctionObject.html>
>  
> <https://www.dealii.org/current/doxygen/deal.II/classVectorFunctionFromScalarFunctionObject.html>
>
> You could pass a lambda function as an argument to the constructor of the 
> function, that uses your FEFieldFunction object, I.e.,
>
> VectorFucntionFromScalarFuncObj my_fun([&](const auto ){ return 
> fe_field_func.value(p);}, 0, 2);
>
> Luca
>
> Il giorno 21 lug 2023, alle ore 22:05, Najwa Alshehri  
> ha scritto:
>
> Dear Daniel, your answer makes sense, Finally it worked.
>
>
> Dear Luca, thank you for your answer. Indeed, the ExactSolution is a class 
> I built myself. I followed your suggestion and it was the magic line that 
> solved the issue. *However,* say that I have in another code instead of 
> using a ExactSolution that I build myself, I used an FeFieldFunction for a 
> reference solution, how would I make this function allocate two components?
>
> Best,
> Najwa
>
> On Friday, July 21, 2023 at 1:54:42 PM UTC+3 luca@gmail.com wrote:
>
>> Dear Najwa, 
>>
>> how many components do the `ExactSolution2` function has? 
>>
>> The error complains about `exact_solution.n_components, n_components` not 
>> being equal, in particular it tells you that `ExactSolution2` only has 
>> one component. 
>>
>> If you built the ExactSolution class yourself, make sure at construction 
>> time you tell Function to allocate two components: 
>>
>> ExactSolution2::ExactSolution2( … ) 
>> : Functions::Function(2) 
>>
>> Best, 
>> Luca. 
>>
>>
>> > On 20 Jul 2023, at 12:15, Najwa Alshehri  wrote: 
>> > 
>> > Hello again, 
>> > 
>> > I have a follow-up question. Does this ComponentSelectFunction work 
>> also with vectors that are not written as blocked vectors? I have applied 
>> it before when I was dealing with blocked vectors and it worked perfectly. 
>> > 
>> > So, I did this, 
>> > const ComponentSelectFunction primal_mask(0,2); 
>> > Later, 
>> > VectorTools::integrate_difference(omega2_dh, 
>> > u_omega2, 
>> > ExactSolution2(), 
>> > cellwise_errors2, 
>> > quadrature, 
>> > VectorTools::L2_norm 
>> > ,_mask); 
>> > const double u2_l2_error = 
>> > VectorTools::compute_global_error(triangulation_omega2, 
>> > cellwise_errors2, 
>> > VectorTools::L2_norm); 
>> > 
>> > And I got the following error!! 
>> > An error occurred in line <455> of file 
>> <../include/deal.II/numerics/vector_tools_integrate_difference.templates.h> 
>> in function 
>> > void dealii::VectorTools::internal::do_integrate_difference(const 
>> dealii::hp::MappingCollection&, const 
>> dealii::DoFHandler&, const InVector&, const 
>> dealii::Function&, OutVector&, 
>> const dealii::hp::QCollection&, const dealii::VectorTools::NormType&, 
>> const dealii::Function*, double) [with int dim = 2; int spacedim 
>> = 2; InVector = dealii::Vector; OutVector = dealii::Vector; 
>> typename InVector::value_type = double] 
>> > The violated condition was: 
>> > 
>> ::dealii::deal_II_exceptions::internals::compare_for_equality(exact_solution.n_components,
>>  
>> n_components) 
>> > Additional information: 
>> > Dimension 1 not equal to 2. 
>> > 
>> > Obviously, we have dimensionality mismatching between u_omega2 and the 
>> exact solution. which means that the mask is not really picking up the 
>> first component of the solution u_omega2. Do you have any suggestions to 
>> fix this issue? 
>> > 
>> > Appreciate your help, 
>> > Najwa 
>> > On Thursday, July 20, 2023 at 11:48:04 AM UTC+3 Najwa Alshehri wrote: 
>> > Thank

Re: [deal.II] Q: Question about extracting part of a vector

2023-07-21 Thread Najwa Alshehri
Dear Daniel, your answer makes sense, Finally it worked.

Dear Luca, thank you for your answer. Indeed, the ExactSolution is a class 
I built myself. I followed your suggestion and it was the magic line that 
solved the issue. *However,* say that I have in another code instead of 
using a ExactSolution that I build myself, I used an FeFieldFunction for a 
reference solution, how would I make this function allocate two components?

Best,
Najwa

On Friday, July 21, 2023 at 1:54:42 PM UTC+3 luca@gmail.com wrote:

> Dear Najwa,
>
> how many components do the `ExactSolution2` function has?
>
> The error complains about `exact_solution.n_components, n_components` not 
> being equal, in particular it tells you that `ExactSolution2` only has 
> one component. 
>
> If you built the ExactSolution class yourself, make sure at construction 
> time you tell Function to allocate two components:
>
> ExactSolution2::ExactSolution2( … ) 
> : Functions::Function(2)
>
> Best,
> Luca.
>
>
> > On 20 Jul 2023, at 12:15, Najwa Alshehri  wrote:
> > 
> > Hello again,
> > 
> > I have a follow-up question. Does this ComponentSelectFunction work also 
> with vectors that are not written as blocked vectors? I have applied it 
> before when I was dealing with blocked vectors and it worked perfectly. 
> > 
> > So, I did this, 
> > const ComponentSelectFunction primal_mask(0,2);
> > Later,
> > VectorTools::integrate_difference(omega2_dh,
> > u_omega2,
> > ExactSolution2(),
> > cellwise_errors2,
> > quadrature,
> > VectorTools::L2_norm
> > ,_mask);
> > const double u2_l2_error =
> > VectorTools::compute_global_error(triangulation_omega2,
> > cellwise_errors2,
> > VectorTools::L2_norm);
> > 
> > And I got the following error!!
> > An error occurred in line <455> of file 
> <../include/deal.II/numerics/vector_tools_integrate_difference.templates.h> 
> in function
> > void dealii::VectorTools::internal::do_integrate_difference(const 
> dealii::hp::MappingCollection&, const 
> dealii::DoFHandler&, const InVector&, const 
> dealii::Function&, OutVector&, 
> const dealii::hp::QCollection&, const dealii::VectorTools::NormType&, 
> const dealii::Function*, double) [with int dim = 2; int spacedim 
> = 2; InVector = dealii::Vector; OutVector = dealii::Vector; 
> typename InVector::value_type = double]
> > The violated condition was: 
> > 
> ::dealii::deal_II_exceptions::internals::compare_for_equality(exact_solution.n_components,
>  
> n_components)
> > Additional information: 
> > Dimension 1 not equal to 2.
> > 
> > Obviously, we have dimensionality mismatching between u_omega2 and the 
> exact solution. which means that the mask is not really picking up the 
> first component of the solution u_omega2. Do you have any suggestions to 
> fix this issue?
> > 
> > Appreciate your help,
> > Najwa
> > On Thursday, July 20, 2023 at 11:48:04 AM UTC+3 Najwa Alshehri wrote:
> > Thank you Daniel for the clear quick answer. I will follow it.
> > 
> > Best,
> > Najwa
> > 
> > On Wednesday, July 19, 2023 at 5:16:14 PM UTC+3 d.arnd...@gmail.com 
> wrote:
> > Najwa,
> > 
> > The documentation of VectorTools::integrate_difference(
> https://www.dealii.org/current/doxygen/deal.II/namespaceVectorTools.html#a676190d2c897ac5da68a9c460fa95832)
>  
> says
> > 
> > 
> > The additional argument weight allows to evaluate weighted norms. The 
> weight function may be scalar, establishing a spatially variable weight in 
> the domain for all components equally. This may be used, for instance, to 
> only integrate over parts of the domain. The weight function may also be 
> vector-valued, with as many components as the finite element: Then, 
> different components get different weights. A typical application is when 
> the error with respect to only one or a subset of the solution variables is 
> to be computed, in which case the other components would have weight values 
> equal to zero. The ComponentSelectFunction class is particularly useful for 
> this purpose as it provides such a "mask" weight. The weight function is 
> expected to be positive, but negative values are not filtered. The default 
> value of this function, a null pointer, is interpreted as "no weighting 
> function", i.e., weight=1 in the whole domain for all vector components 
> uniformly.
> > 
> > Best,
> > Daniel
> > 
> > On Wed, Jul 19, 2023 at 8:01 AM Najwa Alshehri  
> wrote:
> > Dear group members,
> > 
> > I have one question, 
> > 
> > I am trying

Re: [deal.II] Q: Question about extracting part of a vector

2023-07-20 Thread Najwa Alshehri
Hello again,

I have a follow-up question. Does this ComponentSelectFunction work also 
with vectors that are not written as blocked vectors? I have applied it 
before when I was dealing with blocked vectors and it worked perfectly. 

So, I did this, 
const ComponentSelectFunction primal_mask(0,2);
Later,
VectorTools::integrate_difference(omega2_dh,
u_omega2,
ExactSolution2(),
cellwise_errors2,
quadrature,
VectorTools::L2_norm
,_mask);
const double u2_l2_error =
VectorTools::compute_global_error(triangulation_omega2,
cellwise_errors2,
VectorTools::L2_norm);

And I got the following error!!
An error occurred in line <455> of file 
<../include/deal.II/numerics/vector_tools_integrate_difference.templates.h> 
in function
void dealii::VectorTools::internal::do_integrate_difference(const 
dealii::hp::MappingCollection&, const 
dealii::DoFHandler&, const InVector&, const 
dealii::Function&, OutVector&, 
const dealii::hp::QCollection&, const dealii::VectorTools::NormType&, 
const dealii::Function*, double) [with int dim = 2; int spacedim 
= 2; InVector = dealii::Vector; OutVector = dealii::Vector; 
typename InVector::value_type = double]
The violated condition was: 

::dealii::deal_II_exceptions::internals::compare_for_equality(exact_solution.n_components,
 
n_components)
Additional information: 
Dimension 1 not equal to 2.

 Obviously, we have dimensionality mismatching between u_omega2 and the 
exact solution. which means that the mask is not really picking up the 
first component of the solution u_omega2.  Do you have any suggestions to 
fix this issue?

Appreciate your help,
Najwa
On Thursday, July 20, 2023 at 11:48:04 AM UTC+3 Najwa Alshehri wrote:

> Thank you Daniel for the clear quick answer. I will follow it.
>
> Best,
> Najwa
>
> On Wednesday, July 19, 2023 at 5:16:14 PM UTC+3 d.arnd...@gmail.com wrote:
>
>> Najwa,
>>
>> The documentation of VectorTools::integrate_difference(
>> https://www.dealii.org/current/doxygen/deal.II/namespaceVectorTools.html#a676190d2c897ac5da68a9c460fa95832)
>>  
>> says
>>
>>
>> The additional argument weight allows to evaluate weighted norms. The 
>> weight function may be scalar, establishing a spatially variable weight in 
>> the domain for all components equally. This may be used, for instance, to 
>> only integrate over parts of the domain. The weight function may also be 
>> vector-valued, with as many components as the finite element: Then, 
>> different components get different weights. A typical application is when 
>> the error with respect to only one or a subset of the solution variables is 
>> to be computed, in which case the other components would have weight values 
>> equal to zero. The ComponentSelectFunction 
>> <https://www.dealii.org/current/doxygen/deal.II/classComponentSelectFunction.html>
>>  
>> class is particularly useful for this purpose as it provides such a "mask" 
>> weight. The weight function is expected to be positive, but negative values 
>> are not filtered. The default value of this function, a null pointer, is 
>> interpreted as "no weighting function", i.e., weight=1 in the whole domain 
>> for all vector components uniformly.
>>
>> Best,
>> Daniel
>>
>> On Wed, Jul 19, 2023 at 8:01 AM Najwa Alshehri  
>> wrote:
>>
>>> Dear group members,
>>>
>>> I have one question, 
>>>
>>> I am trying to calculate the L2 norm of the error of a solution. In 
>>> particular, I have a *vector of the solution* u_omega and a *FeSystem* 
>>> with two fes.
>>>
>>> I am interested in computing the L2 norm of the error related to the 
>>> first fe using *integrate_difference * function. (Note here u_omega is *not 
>>> a blocked vector*, however, it has two solutions stacked together in 
>>> one vector). Can I extract the solution of the first part from u_omega?
>>>
>>>
>>>
>>> Thank you in advance for your help.
>>>
>>> Best,
>>>
>>> Najwa
>>>
>>> -- 
>>> 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+un...@googlegroups.com.
>>> To view this discussion on the web visit 
>>> https://groups.google.com/d/msgid/dealii/9741338e-2a31-418b-815d-277a5d7cb573n%40googlegroups.co

Re: [deal.II] Q: Question about extracting part of a vector

2023-07-20 Thread Najwa Alshehri
Thank you Daniel for the clear quick answer. I will follow it.

Best,
Najwa

On Wednesday, July 19, 2023 at 5:16:14 PM UTC+3 d.arnd...@gmail.com wrote:

> Najwa,
>
> The documentation of VectorTools::integrate_difference(
> https://www.dealii.org/current/doxygen/deal.II/namespaceVectorTools.html#a676190d2c897ac5da68a9c460fa95832)
>  
> says
>
>
> The additional argument weight allows to evaluate weighted norms. The 
> weight function may be scalar, establishing a spatially variable weight in 
> the domain for all components equally. This may be used, for instance, to 
> only integrate over parts of the domain. The weight function may also be 
> vector-valued, with as many components as the finite element: Then, 
> different components get different weights. A typical application is when 
> the error with respect to only one or a subset of the solution variables is 
> to be computed, in which case the other components would have weight values 
> equal to zero. The ComponentSelectFunction 
> <https://www.dealii.org/current/doxygen/deal.II/classComponentSelectFunction.html>
>  
> class is particularly useful for this purpose as it provides such a "mask" 
> weight. The weight function is expected to be positive, but negative values 
> are not filtered. The default value of this function, a null pointer, is 
> interpreted as "no weighting function", i.e., weight=1 in the whole domain 
> for all vector components uniformly.
>
> Best,
> Daniel
>
> On Wed, Jul 19, 2023 at 8:01 AM Najwa Alshehri  
> wrote:
>
>> Dear group members,
>>
>> I have one question, 
>>
>> I am trying to calculate the L2 norm of the error of a solution. In 
>> particular, I have a *vector of the solution* u_omega and a *FeSystem* 
>> with two fes.
>>
>> I am interested in computing the L2 norm of the error related to the 
>> first fe using *integrate_difference * function. (Note here u_omega is *not 
>> a blocked vector*, however, it has two solutions stacked together in one 
>> vector). Can I extract the solution of the first part from u_omega?
>>
>>
>>
>> Thank you in advance for your help.
>>
>> Best,
>>
>> Najwa
>>
>> -- 
>> 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+un...@googlegroups.com.
>> To view this discussion on the web visit 
>> https://groups.google.com/d/msgid/dealii/9741338e-2a31-418b-815d-277a5d7cb573n%40googlegroups.com
>>  
>> <https://groups.google.com/d/msgid/dealii/9741338e-2a31-418b-815d-277a5d7cb573n%40googlegroups.com?utm_medium=email_source=footer>
>> .
>>
>

-- 
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/3a92adbb-f16d-4ac2-89d7-3ffac4a9bf5fn%40googlegroups.com.


[deal.II] Q: Question about extracting part of a vector

2023-07-19 Thread Najwa Alshehri
Dear group members,

I have one question, 

I am trying to calculate the L2 norm of the error of a solution. In 
particular, I have a *vector of the solution* u_omega and a *FeSystem* with 
two fes.

I am interested in computing the L2 norm of the error related to the first 
fe using *integrate_difference * function. (Note here u_omega is *not a 
blocked vector*, however, it has two solutions stacked together in one 
vector). Can I extract the solution of the first part from u_omega?



Thank you in advance for your help.

Best,

Najwa

-- 
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/9741338e-2a31-418b-815d-277a5d7cb573n%40googlegroups.com.