Re: [deal.II] Solving time dependent heat equation with MPI (PETsc)

2022-01-25 Thread Mark Ma
that's good. feel free to discuss any questions related on this topic, I
will get back to you as soon as possible.

syed ansari  于2022年1月25日周二 14:41写道:

> Thank you very much, Mark, for sharing. This is indeed very helpful for me.
>
> Best Regards,
> Syed Ansari S.
>
> On Tue, Jan 25, 2022 at 5:32 PM Mark Ma  wrote:
>
>> Hi,
>>
>> The code was submitted to the github (
>> https://github.com/dealii/code-gallery/pull/92) but is waiting for
>> checking. Alternatively, you can find the code at my github at:
>>
>> https://github.com/MarkMa1990/code-gallery/tree/master/Distributed_Moving_Laser_Heating
>>
>> Have a nice day,
>> Mark
>>
>> syed ansari  于2022年1月24日周一 09:35写道:
>>
>>> Thank you very much Prof. Bangerth for your kind reply.
>>>
>>> On Sun, Jan 23, 2022 at 10:33 PM Wolfgang Bangerth <
>>> bange...@colostate.edu> wrote:
>>>
>>>> On 1/22/22 2:56 AM, syed ansari wrote:
>>>> >  I tried to find this code in the deal.ii code gallery
>>>> but I am
>>>> > unable to locate it. I request the group members to help me to find
>>>> this code.
>>>>
>>>> I don't think it was ever submitted.
>>>>
>>>> Best
>>>>   W.
>>>>
>>>> --
>>>> 
>>>> Wolfgang Bangerth  email:
>>>> bange...@colostate.edu
>>>> www:
>>>> http://www.math.colostate.edu/~bangerth/
>>>>
>>>> --
>>>> The deal.II project is located at http://www.dealii.org/
>>>> For mailing list/forum options, see
>>>> https://groups.google.com/d/forum/dealii?hl=en
>>>> ---
>>>> You received this message because you are subscribed to the Google
>>>> Groups "deal.II User Group" group.
>>>> To unsubscribe from this group and stop receiving emails from it, send
>>>> an email to dealii+unsubscr...@googlegroups.com.
>>>> To view this discussion on the web visit
>>>> https://groups.google.com/d/msgid/dealii/1e09ca9d-c53d-605b-069d-26d35ea0c25a%40colostate.edu
>>>> .
>>>>
>>> --
>>> 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/5hC7jODg-7k/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/CANaEQub-80n7v-4wDArz_6EyHfzs_8XzxTZvXc3CncNxbBmMWg%40mail.gmail.com
>>> <https://groups.google.com/d/msgid/dealii/CANaEQub-80n7v-4wDArz_6EyHfzs_8XzxTZvXc3CncNxbBmMWg%40mail.gmail.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/CANwa5smXTP_vK5Tb0jBg_zwnC_hBe6PeHkgGLR%2B2KMJX-d-xyA%40mail.gmail.com
>> <https://groups.google.com/d/msgid/dealii/CANwa5smXTP_vK5Tb0jBg_zwnC_hBe6PeHkgGLR%2B2KMJX-d-xyA%40mail.gmail.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 a topic in the
> Google Groups "deal.II User Group" group.
> To unsubscribe from this topic, visit
> https://groups.google.com/d/topic/dealii/5hC7jODg-7k/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/CANaEQuYhXp1KWb4bxinrS8PqGtKSQdgdNO%2BduHXwDBWUw-kJLw%40mail.gmail.com
> <https://groups.google.com/d/msgid/dealii/CANaEQuYhXp1KWb4bxinrS8PqGtKSQdgdNO%2BduHXwDBWUw-kJLw%40mail.gmail.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/CANwa5s%3DX-zG7djZn%3DYY0ZFepCu0wkrDACUjSoMZvjSs-N02KUg%40mail.gmail.com.


Re: [deal.II] Solving time dependent heat equation with MPI (PETsc)

2022-01-25 Thread Mark Ma
Hi,

The code was submitted to the github (
https://github.com/dealii/code-gallery/pull/92) but is waiting for
checking. Alternatively, you can find the code at my github at:
https://github.com/MarkMa1990/code-gallery/tree/master/Distributed_Moving_Laser_Heating

Have a nice day,
Mark

syed ansari  于2022年1月24日周一 09:35写道:

> Thank you very much Prof. Bangerth for your kind reply.
>
> On Sun, Jan 23, 2022 at 10:33 PM Wolfgang Bangerth 
> wrote:
>
>> On 1/22/22 2:56 AM, syed ansari wrote:
>> >  I tried to find this code in the deal.ii code gallery but
>> I am
>> > unable to locate it. I request the group members to help me to find
>> this code.
>>
>> I don't think it was ever submitted.
>>
>> Best
>>   W.
>>
>> --
>> 
>> Wolfgang Bangerth  email: bange...@colostate.edu
>> www: http://www.math.colostate.edu/~bangerth/
>>
>> --
>> The deal.II project is located at http://www.dealii.org/
>> For mailing list/forum options, see
>> https://groups.google.com/d/forum/dealii?hl=en
>> ---
>> You received this message because you are subscribed to the Google Groups
>> "deal.II User Group" group.
>> To unsubscribe from this group and stop receiving emails from it, send an
>> email to dealii+unsubscr...@googlegroups.com.
>> To view this discussion on the web visit
>> https://groups.google.com/d/msgid/dealii/1e09ca9d-c53d-605b-069d-26d35ea0c25a%40colostate.edu
>> .
>>
> --
> 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/5hC7jODg-7k/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/CANaEQub-80n7v-4wDArz_6EyHfzs_8XzxTZvXc3CncNxbBmMWg%40mail.gmail.com
> 
> .
>

-- 
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/CANwa5smXTP_vK5Tb0jBg_zwnC_hBe6PeHkgGLR%2B2KMJX-d-xyA%40mail.gmail.com.


[deal.II] Re: Error during p4est installation

2019-03-11 Thread Mark Ma
Thanks Bruno, it works for my case in a supercomputer.
For some clusters, use module load  before compiling. 

To buid the p4est, edit in p4est-setup.sh to add 
1. module load 
2. export CC=mpicc, CXX=mpicxx 

then use bash ./p4est-setup.sh p4est-x.x.tar.gz install_locations

Mark

在 2015年7月28日星期二 UTC+2下午10:25:41,Bruno Turcksin写道:
>
>
>
> On 07/28/2015 03:06 PM, D. Sarah Stamps wrote: 
> > Hmm, still the same issue. Attached is my config.log file. Can you see 
> > any issues? 
> For some reason CC still uses gcc instead of mpicc. Are you sure there 
> was no typo? Do you do export CC=mpicc and then, in the same terminal 
> run the setup script? What do you get when you type in the terminal: 
>
> echo ${SHELL} 
>
> Bruno 
>

-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.


[deal.II] Compile trilinos with Intel MKL

2018-01-23 Thread Mark Ma
Hello everyone,

For this discussion, I just want to post the working setup in some clusters 
that uses intel MKL instead of BLAS ans LAPACK. Here is the code, 
"thrilino_setup.sh", 

mkdir build
cd build

cmake\
-DTrilinos_ENABLE_Amesos=ON  \
-DTrilinos_ENABLE_Epetra=ON  \
-DTrilinos_ENABLE_Ifpack=ON  \
-DTrilinos_ENABLE_AztecOO=ON \
-DTrilinos_ENABLE_Sacado=ON  \
-DTrilinos_ENABLE_Teuchos=ON \
-DTrilinos_ENABLE_MueLu=ON   \
-DTrilinos_ENABLE_ML=ON  \
-DTrilinos_VERBOSE_CONFIGURE=OFF \
-DTPL_ENABLE_MPI=ON  \





*-DTPL_ENABLE_BLAS:BOOL=ON  \
-DTPL_ENABLE_LAPACK:BOOL=ON  \
-DBLAS_LIBRARY_NAMES="mkl_intel_lp64;mkl_sequential;mkl_core"   
   \
-DBLAS_LIBRARY_DIRS="/opt/software/common/intel/compilers_and_libraries_2017.0.098/linux/mkl/lib/intel64"
 
 \-DLAPACK_LIBRARY_NAMES="" 
 \
-DLAPACK_LIBRARY_DIRS="/opt/software/common/intel/compilers_and_libraries_2017.0.098/linux/mkl/lib/intel64"
 
 \*
-DBUILD_SHARED_LIBS=ON   \
-DCMAKE_VERBOSE_MAKEFILE=OFF \
-DCMAKE_BUILD_TYPE=RELEASE   \
-DCMAKE_INSTALL_PREFIX:PATH=$HOME/installed_prog/trilinos \
../

make install

Best,
Mark

*Reference:*
[1]. https://trilinos.org/pipermail/trilinos-users/2011-May/002502.html

-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.


Re: [deal.II] Solving time dependent heat equation with MPI (PETsc)

2017-10-17 Thread Mark Ma
Sure, no problem.

2017-10-17 15:23 GMT+02:00 Wolfgang Bangerth <bange...@colostate.edu>:

> On 10/17/2017 07:21 AM, Mark Ma wrote:
>
>>
>> PS: I think this could be another tutorial show how to implement this, I
>> would like to contribute if you think it is interesting for Deal.II project.
>>
>
> Yes, I think that would be great. Or, maybe simpler and with less
> documentation requirements: make it a part of the code gallery:
>   http://dealii.org/code-gallery.html
> Want to do that?
>
> Thanks in advance!
>  Wolfgang
>
> --
> 
> Wolfgang Bangerth  email: bange...@colostate.edu
>
>www: http://www.math.colostate.edu/~bangerth/
>
> --
> The deal.II project is located at http://www.dealii.org/
> For mailing list/forum options, see https://groups.google.com/d/fo
> rum/dealii?hl=en
> --- You received this message because you are subscribed to a topic in the
> Google Groups "deal.II User Group" group.
> To unsubscribe from this topic, visit https://groups.google.com/d/to
> pic/dealii/5hC7jODg-7k/unsubscribe.
> To unsubscribe from this group and all its topics, send an email to
> dealii+unsubscr...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.
>

-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.


Re: [deal.II] Solving time dependent heat equation with MPI (PETsc)

2017-10-17 Thread Mark Ma
ep * theta, laplace_matrix_T);


solve_T_run (); // solve and store solution to new_solution_T

// 
old_solution_T = new_solution_T;// used for output, GHOST 
CELLS
old_solution_T_cal = new_solution_T;

output_results (timestep_number);


  }


Finally, the results 

<https://lh3.googleusercontent.com/-2TadguSBFTw/WeYDTK6Do_I/AD4/h7iKKtK4G8o_96qDO84gDU4niUsxGuQZwCLcBGAs/s1600/heat-equation-mpi.gif>

Best,
Mark


PS: I think this could be another tutorial show how to implement this, I 
would like to contribute if you think it is interesting for Deal.II project.

在 2017年10月17日星期二 UTC+2上午2:31:15,Wolfgang Bangerth写道:
>
> On 10/16/2017 01:53 PM, Mark Ma wrote: 
> > 
> > I think this problem lies in the time updating of solution using 
> > old_solution, since the mass_matrix and laplace_matrix have already 
> > eliminated the constraint node, /*mass_matrix_T.vmult 
> > (system_rhs, old_solution_T_cal*//*);*/ is no longer valid for this 
> > case. 
>
> Yes, this sounds reasonable from your description. Have you looked at 
> steps 24, 25, 26 to see how it is done there? I could imagine that you 
> need to think about which degrees of freedom you want eliminated/not 
> eliminated in the matrices you multiply with. To this end, write out the 
> weak form of the problem you need to solve in each time step, and think 
> specifically about the boundary values of the trial and test functions 
> and whether they should be part of the matrix you use on the right hand 
> side or not. It is possible that you need to use a different matrix for 
> the lhs and rhs operations. 
>
> Best 
>   W. 
>
> -- 
>  
> Wolfgang Bangerth  email: bang...@colostate.edu 
>  
> www: http://www.math.colostate.edu/~bangerth/ 
>

-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.


Re: [deal.II] Need for simple applied examples

2017-10-16 Thread Mark Ma
Hi Konstantin,

Very very long time not contact with you, how are you doing? I am sorry I 
really have not too much time that could spend on Deal.II project, the 
Maxwell equation's problem and so on. 

what you propose here I think is right, at least in the view of very 
beginner of deal.II like me, the documentation in deal.ii is not hard to 
start with but really hard to understand and implement more advanced 
techniques like Adaptive mesh, hp refinement, PETsc / Trilinos. Since these 
techniques are introduced in different problems, very beginners may be 
expertise in optics but have to learn these techniques from tutorials of 
fludic-dynamics. Although from mathematics point of view, this is simple. I 
believe if there are some kind of simple examples implementing these 
techniques gradually would greatly decrease the learning and thinking time 
of beginners.

By the way, right now I am calculating temperature in MPI and unfortunately 
stucked. 

Best,
Mark
--
Laboratoire Hubert Curien, UMR CNRS 5516,
Bâtiment F 18 Rue du Professeur Benoît Lauras
42000 Saint-Etienne
FRANCE

在 2015年9月7日星期一 UTC+2下午11:26:55,Konstantin Ladutenko写道:
>
> Hi,
>
> Sorry, it was hard for me to find enough time to get really started with 
> dell. Few things I was able to learn this summer you can find in this group 
> and in updated step-6.
>
> The next possible step for me is to solve Mie problem using deal.ii and 
> compare against analytic solution. It something similar to H. Fahs 
> self-written code you can find at http://hfahs.free.fr/software.html (and 
> related bunch of publications http://hfahs.free.fr/pubs.html ). For sure 
> I will present it here as soon as I have something to solve the problem of 
> EM scattering correctly (and probably after that it can be added some 
> nonlinearity for plasmonic nanoparticles). 
>
> I also found a book of Peter Monk "FEM for Maxwell equations" that is 
> rather hard for me to read, however, I hope, it will give me answers to all 
> my basic questions. There is also "Introduction to the Finite Element 
> Method in Electromagnetics" by
> Polycarpou which is much easier to read from my "applied" point of view.
>
> Best regards,
> Kostya
>
>
> пн, 7 сент. 2015 г. в 21:21, JT Hu :
>
>> Hi Kostya,
>>
>> It has been a year since you posted this.  Do you have any updates on 
>> developing Deal II methods for plasmonics/photonics calculations to share? 
>> If you have published any related work discuss the methods and coding, I 
>> will also be glad to read and cite!
>>
>> Best,
>> Jingtian
>>
>>
>> On Tuesday, September 16, 2014 at 4:44:42 AM UTC-5, Konstantin Ladutenko 
>> wrote:
>>
>>> Hi Wolfgang,
>>>
>>>  Thank you for your reply! Deal.II documentation is great, I really like 
>>> your lectures, I will recommend them to my students.
>>>
>>> 1)  My goal is to correctly solve my problem as fast as possible, time 
>>> needed to learn new tool included.
>>>
>>> 2) My choice is open source. We have few Comsol licenses usually not 
>>> available for students (taken out for 'great and terrible' top science 
>>> projects), this also leads to some troubles than you apply to use external 
>>> cluster for computations. There is a number of other open source FEM codes 
>>> (Elmer, GetFEM++, Moose, etc.), however as to me deal.ii seems to be a 
>>> mature tool which I am learning with hope to use it many years after that.
>>>
>>> 3) I suppose that almost any project can benefit from large user base. 
>>> It seems that deal.ii can solve many problems in different research areas. 
>>> However average researcher with need of calculations will prefer to use 
>>> matlab (or octave, scilab), may be some python codes (e.g. with help of 
>>> sage) and sometimes some ancient Fortran codes - too large and awful to be 
>>> maintained. GUI programs (Comsol, Ansys etc) are also an option.
>>>
>>> It is a good way to start from basic elements of deal.ii, as it is done 
>>> in stepXX tutorials and to go to some real world problems step by step. 
>>> However you will need few weeks before your can start with your problem. 
>>> Many possible users of deal.ii do not have so much time to be spend without 
>>> any practical results. May be I am wrong, but a number of simple applied 
>>> examples can be used as a kick start entry point to deal.ii. You can use 
>>> them to copy-paste as building bricks for your exact (may be simplified) 
>>> problem. It will not be solved in the best way deal.ii can provide,  but it 
>>> will motivate to learn more details about how it should be done in deal.ii 
>>> way.
>>>
>>> 4) My current research area of interest is metamaterials (with nonlinear 
>>> effects coming soon), so Maxwell equations are one the first place. I am 
>>> mostly interested in transient problems, however after 20 video lectures I 
>>> am still not sure how i should treat them in deal.ii). This is my entry 
>>> list of examples for electromagnetism (all of them can be 

Re: [deal.II] Solving time dependent heat equation with MPI (PETsc)

2017-10-16 Thread Mark Ma
source code are attached.

在 2017年10月16日星期一 UTC+2下午9:14:13,Wolfgang Bangerth写道:
>
> On 10/16/2017 09:23 AM, Mark Ma wrote: 
> > 
> > It almost looks to me like you're applying both Dirichlet values (via 
> > the ConstraintMatrix) and Neumann boundary values (via a boundary 
> > integral). But then I haven't taken a look at the code, so I can't 
> > really say for sure. 
> > 
> > 
> > I think I only applied Dirichlet values via ConstraintMatrix. There 
> > is no boundary integrals, so I believe Neumann BC is not 
> > implemented. 
>
> My guess came from the fact that a Neumann boundary condition would act 
> like a boundary heat source. It looks like you have a heat source in the 
> last layer of cells around the boundary. 
>
> I don't know whether that is true, but it's a place I'd try to 
> investigate in debugging. 
>
>
> > But may be you are right since for most of the case, heat equations 
> > are solved assuming an adiabatic process, that is dT/dx = 0 at 
> > boundaries. This zero values is automatically satisfied. 
>
> That's a different question. Whether one may *physically* want a 
> different equation is besides the point -- you have chosen one equation 
> you'd like to solve, but the numerical solution is wrong. It's 
> worthwhile trying to understand why that is before thinking about 
> *other* equations. 
>
>
> > I think the problem comes from the way of setting Boundary values 
> > using constraint method, which in dealii it makes some kind of 
> > constraints at vertex of boundaries, i.e. x22 = x1/2 + x2/2, 
>
> Hm, that would be the constraint of a hanging node. Do you have hanging 
> nodes in your problem? If you do, may I suggest to simplify the problem 
> and see if you have the same problem on uniform meshes? 
>
> Best 
>   W. 
>
>
> -- 
>  
> Wolfgang Bangerth  email: bang...@colostate.edu 
>  
> www: http://www.math.colostate.edu/~bangerth/ 
>

-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.
/* -
 *
 * Copyright (C) 1999 - 2016 by the deal.II authors
 *
 * This file is part of the deal.II library.
 *
 * The deal.II library is free software; you can use it, redistribute
 * it, and/or modify it under the terms of the GNU Lesser General
 * Public License as published by the Free Software Foundation; either
 * version 2.1 of the License, or (at your option) any later version.
 * The full text of the license can be found in the file LICENSE at
 * the top level of the deal.II distribution.
 *
 * -

 *
 * Author: Wolfgang Bangerth, University of Heidelberg, 1999
 */


// @sect3{Include files}

// The first few (many?) include files have already been used in the previous
// example, so we will not explain their meaning here again.
//#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
//#include 
#include 
//#include 
//#include 

/*
#include 
#include 
#include 
#include 
*/

#include 
#include 
#include 
#include 


#include 
#include 

#include 
#include 
#include 
#include 

#include 
#include 
#include 

// This is new, however: in the previous example we got some unwanted output
// from the linear solvers. If we want to suppress it, we have to include this
// file and add a single line somewhere to the program (see the main()
// function below for that):
#include 

#include 

#include 

// The final step, as in previous programs, is to import all the deal.II class
// and function names into the global namespace:
using namespace dealii;


// @sect3{The Step4 class template}

// This is again the same Step4 class as in the previous
// example. The only difference is that we have now declared it as a class
// with a template parameter, and the template parameter is of course the
// spatial dimension in which we would like to solve the Laplace equation. Of
// course, several of the member variables depend on this dimension as well,
// in particular the Triangulation class, which has to represent
// quadrilaterals or hexahedra, respectively. Apart from this, everything is
// as before.
template 
class Step4
{
public:
  Step4 ();
  void run ();

private:
  voi

Re: [deal.II] Solving time dependent heat equation with MPI (PETsc)

2017-10-16 Thread Mark Ma


> On 10/16/2017 09:23 AM, Mark Ma wrote: 
> > 
> > It almost looks to me like you're applying both Dirichlet values (via 
> > the ConstraintMatrix) and Neumann boundary values (via a boundary 
> > integral). But then I haven't taken a look at the code, so I can't 
> > really say for sure. 
> > 
> > 
> > I think I only applied Dirichlet values via ConstraintMatrix. There 
> > is no boundary integrals, so I believe Neumann BC is not 
> > implemented. 
>
> My guess came from the fact that a Neumann boundary condition would act 
> like a boundary heat source. It looks like you have a heat source in the 
> last layer of cells around the boundary. 
>
> I don't know whether that is true, but it's a place I'd try to 
> investigate in debugging. 
>
>
Yes, it looks like so. But then I comment the time updating code of rhs 
source by the following code, this boundary heat still appears(see figs 
below)
 //
  //time dependent

//assign right hand side
mass_matrix_T.vmult (system_rhs, old_solution_T_cal);
laplace_matrix_T.vmult (tmp, old_solution_T_cal);
system_rhs.add (-time_step * (1-theta), tmp);
/*
assemble_rhs_T (time);
forcing_terms = dynamic_rhs_T;
forcing_terms *= time_step * theta;

assemble_rhs_T (time - time_step); 
tmp = dynamic_rhs_T;
forcing_terms.add (time_step*(1-theta),tmp);
system_rhs.add (1,forcing_terms);

*/

//assign system matrix
system_matrix.copy_from (mass_matrix_T);
system_matrix.add (time_step * theta, laplace_matrix_T);

I think this problem lies in the time updating of solution using 
old_solution, since the mass_matrix and laplace_matrix have already 
eliminated the constraint node, 
*mass_matrix_T.vmult (system_rhs, old_solution_T_cal**);* is no longer 
valid for this case. This may result into ZEROs at all constrained nodes, 
so at boundaries these zero values then propagate into center, as can be 
seen from the figure below. For time updating of solution using 
old_solution, maybe tutorial-32 could give some ideas. However to 
completely understand step-32 may spend some more time for me. 

<https://lh3.googleusercontent.com/-yfvZmZDrpHE/WeUMm6TyO-I/ADo/b0oiTCd5_HMGIsAe9LlZ-RPYF5nWu9ygwCLcBGAs/s1600/tem_no_source.gif>




 

>
> > But may be you are right since for most of the case, heat equations 
> > are solved assuming an adiabatic process, that is dT/dx = 0 at 
> > boundaries. This zero values is automatically satisfied. 
>
> That's a different question. Whether one may *physically* want a 
> different equation is besides the point -- you have chosen one equation 
> you'd like to solve, but the numerical solution is wrong. It's 
> worthwhile trying to understand why that is before thinking about 
> *other* equations. 
>
>
> > I think the problem comes from the way of setting Boundary values 
> > using constraint method, which in dealii it makes some kind of 
> > constraints at vertex of boundaries, i.e. x22 = x1/2 + x2/2, 
>
> Hm, that would be the constraint of a hanging node. Do you have hanging 
> nodes in your problem? If you do, may I suggest to simplify the problem 
> and see if you have the same problem on uniform meshes? 
>
>
there is no handing node in my code since I used a uniformed mesh. And also 
I tried this, results show no difference.

I am really appreciate your patience to look for this question during your 
busy schedule.

Thanks and have a nice day,
Mark
-
Laboratoire Hubert Curien, UMR CNRS 5516,
Bâtiment F 18 Rue du Professeur Benoît Lauras
42000 Saint-Etienne
FRANCE

-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.


Re: [deal.II] Solving time dependent heat equation with MPI (PETsc)

2017-10-16 Thread Mark Ma

在 2017年10月16日星期一 UTC+2下午4:31:25,Wolfgang Bangerth写道:
>
>
> > So when you visualize the solution, the error is at the boundary but 
> it looks 
> > correct in the interior? 
> > 
> > Yes, it is. Thereafter, the error at boundaries does propagate and then 
> > interfere with the interior. 
>
>  From the pictures you posted in a follow-up, it looks like the boundary 
> values are actually correct, but that the problem starts one layer of 
> cells 
> into the domain. Do I see this correctly? 
>
>
Yes, exactly. And if I use the BC values as 0, then it seems OK now (see 
figs below). What if I want to apply some value other than 0, this problem 
seems annoying.
 

>
> >  > Is there a simple way to update the RHS of old value using 
> something 
> > simple 
> >  > like vmult? 
> > 
> > What is the equation you are trying to implement? 
> > 
> > I am trying just the Heat equation, 
> > rho * C * dT/dt = nabla *(kappba nalbaT) + f(x,y,t); 
> > Boundary condition is, 
> > T(x,y,z,t) = 300K, at boundaries. 
> I meant the "update the RHS of old value" equation. 
>   
>

sorry, I mistook your idea. The RHS function is,
f(x,y,z,t) = I0 * exp(-((x-v*t)^2 + (y-v*t)^2 +(z-v*t)^2)/2c^2), which 
stands for absorption of the laser.





Boundary value = 0, initial values = 0


-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.


[deal.II] Re: Solving time dependent heat equation with MPI (PETsc)

2017-10-16 Thread Mark Ma



<https://lh3.googleusercontent.com/-6nIRhm7FQYU/WeR-GMJGCuI/ADA/HbpSoiEtK5cHx02THJj9izq1lQE2jUlvgCLcBGAs/s1600/temperature.gif>


<https://lh3.googleusercontent.com/-6nIRhm7FQYU/WeR-GMJGCuI/ADA/HbpSoiEtK5cHx02THJj9izq1lQE2jUlvgCLcBGAs/s1600/temperature.gif>


<https://lh3.googleusercontent.com/-RwHcXcK15qY/WeR-JkzMwNI/ADE/xHC76YIouRAWIAZmoJ3kkUMZpZBmi0rCACLcBGAs/s1600/Temperature..png>

Projection of initial value, time 000
<https://lh3.googleusercontent.com/-RwHcXcK15qY/WeR-JkzMwNI/ADE/xHC76YIouRAWIAZmoJ3kkUMZpZBmi0rCACLcBGAs/s1600/Temperature..png>


<https://lh3.googleusercontent.com/-RwHcXcK15qY/WeR-JkzMwNI/ADE/xHC76YIouRAWIAZmoJ3kkUMZpZBmi0rCACLcBGAs/s1600/Temperature..png>


<https://lh3.googleusercontent.com/-t9fNhCEb15E/WeR-MXkwLdI/ADI/T11IwROMgWswp4z8wp3O7urkAL9F5BgQQCLcBGAs/s1600/Temperature.0001.png>

time 001
<https://lh3.googleusercontent.com/-t9fNhCEb15E/WeR-MXkwLdI/ADI/T11IwROMgWswp4z8wp3O7urkAL9F5BgQQCLcBGAs/s1600/Temperature.0001.png>


<https://lh3.googleusercontent.com/-gYujhXn8yJQ/WeR-PY7fXXI/ADM/dL2wiA4zUuEjkvDUzuIwKwjCOoyMebFKACLcBGAs/s1600/Temperature.0002.png>

time 002
<https://lh3.googleusercontent.com/-gYujhXn8yJQ/WeR-PY7fXXI/ADM/dL2wiA4zUuEjkvDUzuIwKwjCOoyMebFKACLcBGAs/s1600/Temperature.0002.png>


<https://lh3.googleusercontent.com/-t9fNhCEb15E/WeR-MXkwLdI/ADI/T11IwROMgWswp4z8wp3O7urkAL9F5BgQQCLcBGAs/s1600/Temperature.0001.png>


<https://lh3.googleusercontent.com/-RwHcXcK15qY/WeR-JkzMwNI/ADE/xHC76YIouRAWIAZmoJ3kkUMZpZBmi0rCACLcBGAs/s1600/Temperature..png>



在 2017年10月10日星期二 UTC+2下午4:40:34,Mark Ma写道:
>
>
> Dear experienced deal.II users and developers, 
>
> I want to solve a heat equation in the time domain with distributed memory 
> using MPI, but the results are incorrect. In order to do so, I reference 
> tutorial step-23 for time updating method and step-40 for implementing MPI. 
> May I ask whether my boundary condition is right or not? Should we do 
> compress() after apply_boundary_values()? Thanks in advance!
>
> With best regards,
> Mark
>
>
>
>
>
> In PETsc or Trilinos, how could we set the boundary conditions and initial 
> conditions?
>
>
> Since old_solution_T would be read in other places (not shown here), it is 
> initialized with ghost cells; while old_solutuon_T_cal is only used for 
> written, it does not have ghost cells. see code as following,
> //
> old_solution_T_cal.reinit (locally_owned_dofs,mpi_communicator);
> old_solution_T.reinit 
> (locally_owned_dofs,locally_relevant_dofs,mpi_communicator);
> //
> ..
>
>
>   template 
>   void ThermalDiffusion::run ()
>   {
>
>
> setup_system();
> assemble_system();
>
> 
> VectorTools::project (dof_handler, constraints, QGauss(degree),
>   InitialValues_T(),
>   old_solution_T_cal);
>
> old_solution_T = old_solution_T_cal;
>
> LA::MPI::Vector tmp (locally_owned_dofs,mpi_communicator);
> LA::MPI::Vector forcing_terms (locally_owned_dofs,mpi_communicator);
>
> for (timestep_number=1, time=time_step;
>  time<= global_simulation_end_time;
>  time+=time_step, ++timestep_number)
>   {
> pcout << "Time step " << timestep_number
>   << " at t=" << time
>   << std::endl;
>
>
> //---
> //run to solve T
> //---
>
>   //
>   //time dependent
>
> //assign right hand side
> mass_matrix_T.vmult (system_rhs, old_solution_T_cal);
>
> laplace_matrix_T.vmult (tmp, old_solution_T_cal);
> system_rhs.add (-time_step * (1-theta), tmp);
>
> assemble_rhs_T (time);
> forcing_terms = dynamic_rhs_T;
> forcing_terms *= time_step * theta;
>
> assemble_rhs_T (time - time_step); 
> tmp = dynamic_rhs_T;
> forcing_terms.add (time_step*(1-theta),tmp);
> system_rhs.add (1,forcing_terms);
>
> //assign system matrix
> system_matrix.copy_from (mass_matrix_T);
> system_matrix.add (time_step * theta, laplace_matrix_T);
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
> * {  BoundaryValues_Temperature 
> boundary_values_function_T;  //boundary_values_function.set_time 
> (time);  std::map<types::global_dof_index,double> 
> boundary_values_T;  VectorTools::interpolate_boundary_values 
> (dof_handler, 

Re: [deal.II] Solving time dependent heat equation with MPI (PETsc)

2017-10-16 Thread Mark Ma


> On 10/15/2017 12:41 PM, Mark Ma wrote: 
> > 
> > Now the projection of initial values (rewrite the code by manually 
> assemble 
> > the matrix and system_rhs and calculate) run OK, but the time updating 
> of T is 
> > not correct, same phenomenon appears. I believe this may arise from the 
> fact 
> > that direct using matrix vmult (i.e. 
> > | 
> > mass_matrix_T.vmult (system_rhs,old_solution_T_cal); 
> > | 
> > 
> > ) instead of assembling and distribute_local_to_global again may ignore 
> > eliminating the constraint points in matrix or vector, when using 
> > constriantMatrix and interporate_boundary_values to apply the boundary 
> > condition, I am now checking of this. 
>
> So when you visualize the solution, the error is at the boundary but it 
> looks 
> correct in the interior? 
>
> Yes, it is. Thereafter, the error at boundaries does propagate and then 
interfere with the interior. 

>
> > Is there a simple way to update the RHS of old value using something 
> simple 
> > like vmult? 
>
> What is the equation you are trying to implement? 
>
> I am trying just the Heat equation,
rho * C * dT/dt = nabla *(kappba nalbaT) + f(x,y,t);
Boundary condition is,
T(x,y,z,t) = 300K, at boundaries.
 

> Best 
>   W. 
>
> -- 
>  
> Wolfgang Bangerth  email: bang...@colostate.edu 
>  
> www: http://www.math.colostate.edu/~bangerth/ 
>
>

-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.


Re: [deal.II] Solving time dependent heat equation with MPI (PETsc)

2017-10-15 Thread Mark Ma
Prof. W. Bangerth,

Please find attached source code, where the projection part is ok now, but 
time updating of solution part results is incorrect I believe (see boundary 
results). This confused me a few weeks and still did not find a solution. 

Best,
Mark

在 2017年10月13日星期五 UTC+2下午2:55:04,Wolfgang Bangerth写道:
>
> On 10/13/2017 02:06 AM, Mark Ma wrote: 
> > 
> > later, I changed the control into 
> > | 
> > | 
> > SolverControlsolver_control 
> (5*system_rhs.size(),1e-12*system_rhs.l2_norm()) 
> > | 
> > | 
> > this works well for structure size of um or nm. I think previous setting 
> may 
> > lead to a loss of precision so that the results are always incorrect. 
>
> Yes, indeed -- using a tolerance relative to the size of the rhs vector is 
> important. 
>
> Best 
>   W. 
>
> -- 
>  
> Wolfgang Bangerth  email: bang...@colostate.edu 
>  
> www: http://www.math.colostate.edu/~bangerth/ 
>
>

-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.
/* -
 *
 * Copyright (C) 1999 - 2016 by the deal.II authors
 *
 * This file is part of the deal.II library.
 *
 * The deal.II library is free software; you can use it, redistribute
 * it, and/or modify it under the terms of the GNU Lesser General
 * Public License as published by the Free Software Foundation; either
 * version 2.1 of the License, or (at your option) any later version.
 * The full text of the license can be found in the file LICENSE at
 * the top level of the deal.II distribution.
 *
 * -

 *
 * Author: Wolfgang Bangerth, University of Heidelberg, 1999
 */


// @sect3{Include files}

// The first few (many?) include files have already been used in the previous
// example, so we will not explain their meaning here again.
//#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
//#include 
#include 
//#include 
//#include 

/*
#include 
#include 
#include 
#include 
*/

#include 
#include 
#include 
#include 


#include 
#include 

#include 
#include 
#include 
#include 

#include 
#include 
#include 

// This is new, however: in the previous example we got some unwanted output
// from the linear solvers. If we want to suppress it, we have to include this
// file and add a single line somewhere to the program (see the main()
// function below for that):
#include 

#include 

#include 

// The final step, as in previous programs, is to import all the deal.II class
// and function names into the global namespace:
using namespace dealii;


// @sect3{The Step4 class template}

// This is again the same Step4 class as in the previous
// example. The only difference is that we have now declared it as a class
// with a template parameter, and the template parameter is of course the
// spatial dimension in which we would like to solve the Laplace equation. Of
// course, several of the member variables depend on this dimension as well,
// in particular the Triangulation class, which has to represent
// quadrilaterals or hexahedra, respectively. Apart from this, everything is
// as before.
template 
class Step4
{
public:
  Step4 ();
  void run ();

private:
  void make_grid ();
  void setup_system();
  void assemble_system ();

  void assemble_rhs_T(double rhs_time);
  void solve_T ();
  void solve_T_run ();
  void output_results (int output_num) const;
  
  MPI_Commmpi_communicator;


  parallel::distributed::Triangulation   triangulation;
  FE_Qfe;
  DoFHandler  dof_handler;

  ConstraintMatrix constraints;

//  SparsityPattern  sparsity_pattern;

  // system_matrix firstly used for projection of
  //init_T
  TrilinosWrappers::SparseMatrix system_matrix_T;
  TrilinosWrappers::SparseMatrix mass_matrix_T;
  TrilinosWrappers::SparseMatrix laplace_matrix_T;

  // with ghost cells
  //old_solution_T
  TrilinosWrappers::MPI::Vector   old_solution_T; 

  // only locally owned cells
  //old_solution_T_cal
  //system_rhs
  TrilinosWrappers::MPI::Vector   old_solution_T_cal;
  TrilinosWrappers::MPI::Vector   new_solution_T;
  TrilinosWrappers::MPI::Vector   system_rhs_T;

// also store right ha

[deal.II] Re: Error when applying initial values to MPI::Vector in multiple dimensions

2017-10-15 Thread Mark Ma
Dear Maxi,

For projecting initial values in MPI (project function may fail using PETsc 
or Trilinos), it is more convenient by direct solving the equation:
Mass_matrix*Solution =(InitialvaluesFunc,phi)
this part in MPI is not hard, you could do it very quickly I think.

Best,
Mark


在 2017年10月13日星期五 UTC+2上午11:05:01,Maxi Miller写道:
>
> I try to apply initial values to a vector defined as 
> *LinearAlgebraTrilinos::MPI::Vector 
> *using 
> VectorTools::project (dof_handler, hanging_node_constraints,
>  QGauss(fe.degree+1),
>  InitialValues(),
>  local_solution);
>
>
>
> When initializing the variable fe (as FESystem) with one or two 
> components, it works fine. For more than two components I get the error
>
>  
> An error occurred in line <1366> of file 
> <~/Downloads/dealii/include/deal.II/numerics/vector_tools.templates.h> in 
> function 
> void dealii::VectorTools::{anonymous}::project(const 
> dealii::Mapping&, const dealii::DoFHandler&, const 
> dealii::ConstraintMatrix&, const dealii::Quadrature&, const 
> dealii::Function&, VectorType&, bool, 
> const dealii
> ::Quadrature<(dim - 1)>&, bool) [with VectorType = 
> dealii::TrilinosWrappers::MPI::Vector; int dim = 2; typename 
> VectorType::value_type = double] 
> The violated condition was:  
> (dynamic_cast* > 
> (&(dof.get_triangulation()))==nullptr) 
> Additional information:  
> You are trying to use functionality in deal.II that is currently not 
> implemented. In many cases, this indicates that there simply didn't appear 
> much of a need for it, or that the author of the original code did not have 
> the time to implement a particular case. If you
>  hit this exception, it is therefore worth the time to look into the code to 
> find out whether you may be able to implement the missing functionality. If 
> you do, please consider providing a patch to the deal.II development sources 
> (see the deal.II website on how to contri
> bute). 
>  
> Stacktrace: 
> --- 
> #0  /opt/dealII/lib/libdeal_II.g.so.9.0.0-pre:  
> #1  /opt/dealII/lib/libdeal_II.g.so.9.0.0-pre: void 
> dealii::VectorTools::project<2, dealii::TrilinosWrappers::MPI::Vector, 
> 2>(dealii::DoFHandler<2, 2> const&, dealii::ConstraintMatrix const&, 
> dealii::Quadrature<2> const&, dealii::Function<2, 
> dealii::TrilinosWrappers::MPI
> ::Vector::value_type> const&, dealii::TrilinosWrappers::MPI::Vector&, bool, 
> dealii::Quadrature<(2)-(1)> const&, bool) 
> #2  ./main: Step15::MinimalSurfaceProblem<2>::run() 
> #3  ./main: main 
>  
>  
> [linux-lb8c:15830] *** Process received signal *** 
> [linux-lb8c:15830] Signal: Aborted (6) 
> [linux-lb8c:15830] Signal code:  (-6) 
> [linux-lb8c:15830] [ 0] /lib64/libpthread.so.0(+0x12270)[0x7f294a477270] 
> [linux-lb8c:15830] [ 1] /lib64/libc.so.6(gsignal+0x110)[0x7f2946c1f0d0] 
> [linux-lb8c:15830] [ 2] /lib64/libc.so.6(abort+0x151)[0x7f2946c206b1] 
> [linux-lb8c:15830] [ 3] 
> /opt/dealII/lib/libdeal_II.g.so.9.0.0-pre(+0x6b9e5d1)[0x7f295b49e5d1] 
> [linux-lb8c:15830] [ 4] 
> /opt/dealII/lib/libdeal_II.g.so.9.0.0-pre(_ZN6dealii18deal_II_exceptions9internals5abortERKNS_13ExceptionBaseE+0x1a)[0x7f295b49edaf]
>  
> [linux-lb8c:15830] [ 5] 
> /opt/dealII/lib/libdeal_II.g.so.9.0.0-pre(_ZN6dealii18deal_II_exceptions9internals11issue_errorINS_18StandardExceptions17ExcNotImplementedEEEvNS1_17ExceptionHandlingEPKciS7_S7_S7_T_+0x98)[0x7f2957373ea1]
>  
> [linux-lb8c:15830] [ 6] 
> /opt/dealII/lib/libdeal_II.g.so.9.0.0-pre(+0x3f38e23)[0x7f2958838e23] 
> [linux-lb8c:15830] [ 7] 
> /opt/dealII/lib/libdeal_II.g.so.9.0.0-pre(_ZN6dealii11VectorTools7projectILi2ENS_16TrilinosWrappers3MPI6VectorELi2EEEvRKNS_10DoFHandlerIXT_EXT1_EEERKNS_16ConstraintMatrixERKNS_10QuadratureIXT_EEERKNS_8FunctionIXT1_ENT0_10value_typeEEERSH_bRKNSC_IX
> miT_Li1b+0x2f)[0x7f295894906e] 
> [linux-lb8c:15830] [ 8] 
> ./main(_ZN6Step1521MinimalSurfaceProblemILi2EE3runEv+0xc08)[0x420d08] 
> [linux-lb8c:15830] [ 9] ./main(main+0x3c)[0x414ad0] 
> [linux-lb8c:15830] [10] 
> /lib64/libc.so.6(__libc_start_main+0xea)[0x7f2946c09f4a] 
> [linux-lb8c:15830] [11] ./main(_start+0x2a)[0x41477a] 
> [linux-lb8c:15830] *** End of error message *** 
> Abgebrochen (Speicherabzug geschrieben)
>
> when running in debug mode. It runs fine in release mode. Why does that 
> happen for more than two components, and how can I fix/circumvent that? Or 
> did I (again) forget something? 
>
> My minimal example is attached, the behaviour happens when setting 
> NUM_COMPONENTS via 
>
> #define NUM_COMPONENTS 100
>
> to a value larger than 2.
>
>
> 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 

Re: [deal.II] Solving time dependent heat equation with MPI (PETsc)

2017-10-15 Thread Mark Ma
Prof. Wolfgang Bangerth,

Now the projection of initial values (rewrite the code by manually assemble 
the matrix and system_rhs and calculate) run OK, but the time updating of T 
is not correct, same phenomenon appears. I believe this may arise from the 
fact that direct using matrix vmult (i.e. 
mass_matrix_T.vmult (system_rhs, old_solution_T_cal);

) instead of assembling and distribute_local_to_global again may ignore 
eliminating the constraint points in matrix or vector, when using 
constriantMatrix and interporate_boundary_values to apply the boundary 
condition, I am now checking of this. 

Is there a simple way to update the RHS of old value using something simple 
like vmult?

Best,
Mark

//---
//Updating of T in time domain--the old way, working in non-mpi version
//---

  //
  //time dependent

//assign right hand side
mass_matrix_T.vmult (system_rhs, old_solution_T_cal);

laplace_matrix_T.vmult (tmp, old_solution_T_cal);
system_rhs.add (-time_step * (1-theta), tmp);

assemble_rhs_T (time);
forcing_terms = dynamic_rhs_T;
forcing_terms *= time_step * theta;

assemble_rhs_T (time - time_step); 
tmp = dynamic_rhs_T;
forcing_terms.add (time_step*(1-theta),tmp);
system_rhs.add (1,forcing_terms);

//assign system matrix
system_matrix.copy_from (mass_matrix_T);
system_matrix.add (time_step * theta, laplace_matrix_T);




在 2017年10月13日星期五 UTC+2下午4:43:30,Wolfgang Bangerth写道:
>
> On 10/13/2017 08:39 AM, Lucas Campos wrote: 
> > 
> > In general, using MatrixTools::apply_boundary_values() is not the 
> way to go 
> > with MPI programs. Rather, use a ConstraintMatrix and incorporate 
> the 
> > boundary 
> > values into the same object as you do with hanging node constraints. 
> > 
> > 
> > This is the way to go due to correctness, or in the sense of 
> scalability? 
>
> MatrixTools::apply_boundary_values() needs access to the elements of the 
> matrix because it wants to modify elements after they have already been 
> written into the matrix. That is already difficult if the matrix is owned 
> by 
> PETSc or Trilinos -- we can get access to these elements, but it is not 
> efficient to do so. 
>
> But the bigger issue is that the function wants to access elements not 
> only 
> for the rows of constrained DoFs, but also for the columns. That means 
> that 
> you may have to access elements that are actually stored on other 
> processors 
> -- something that can not be done efficiently. Consequently, 
> MatrixTools::apply_boundary_values() does not attempt to eliminate columns 
> of 
> the matrix, and you will end up with a non-symmetric matrix even if your 
> problem is symmetric. 
>
> It is better to use the approach via ConstraintMatrix that deals with 
> entries 
> before they even get into the matrix (and therefore in particular before 
> matrix entries are sent to other processors). 
>
> Best 
>   W. 
>
> -- 
>  
> Wolfgang Bangerth  email: bang...@colostate.edu 
>  
> www: http://www.math.colostate.edu/~bangerth/ 
>
>

-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.


Re: [deal.II] Solving time dependent heat equation with MPI (PETsc)

2017-10-13 Thread Mark Ma
It is actually my stupid mistake in solvercontrol,

I used
SolverControl solver_control (dof_handler.n_dofs(),1e-12)
This works when the geometry size in the order of 1, but fails at 1e-6, or
even 1e-9. Here I actually want to make a micrometer or nanometer size
structure.

later, I changed the control into
SolverControl solver_control (5*system_rhs.size(),1e-12*system_rhs.l2_norm()
)
this works well for structure size of um or nm. I think previous setting
may lead to a loss of precision so that the results are always incorrect.

Thanks for your ideas and have a nice day,
Mark

2017-10-12 14:29 GMT+02:00 Wolfgang Bangerth <bange...@colostate.edu>:

> On 10/12/2017 03:02 AM, Mark Ma wrote:
>
>>
>> Thanks for your useful advice. Finally I solved this problem.
>>
>
> So what was the problem? Maybe we can learn from your example?
>
>
> Best
>  W.
>
> --
> 
> Wolfgang Bangerth  email: bange...@colostate.edu
>www: http://www.math.colostate.edu/~bangerth/
>
> --
> The deal.II project is located at http://www.dealii.org/
> For mailing list/forum options, see https://groups.google.com/d/fo
> rum/dealii?hl=en
> --- You received this message because you are subscribed to a topic in the
> Google Groups "deal.II User Group" group.
> To unsubscribe from this topic, visit https://groups.google.com/d/to
> pic/dealii/5hC7jODg-7k/unsubscribe.
> To unsubscribe from this group and all its topics, send an email to
> dealii+unsubscr...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.
>

-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.


Re: [deal.II] Solving time dependent heat equation with MPI (PETsc)

2017-10-12 Thread Mark Ma
Hi Prof. Wolfgang Bangerth,

Thanks for your useful advice. Finally I solved this problem.

Best,
Mark

2017-10-10 19:50 GMT+02:00 Wolfgang Bangerth :

> On 10/10/2017 11:40 AM, Mark wrote:
>
>>
>> I really appreciate for your quick reply. Now I do as what you said,
>> using VectorTools in setup_system(). It seems working, but there is still
>> one question that I always find very large values at points near
>> boundaries, see attached figures. Is there a way to solve this?
>>
>
> There are so many possible reasons for this that I can't even speculate.
> Play with what you have to find out what the common theme is -- how do
> things change if you try different boundary conditions, initial conditions,
> right hand sides, time steps, etc. What affects the problem and what
> doesn't? Can you make the problem simpler (e.g., use zero boundary values
> and zero initial conditions)? Does this only happen in parallel or also if
> you run the program on a single processor?
>
> This is sort of the mental list that I go through when debugging these
> kinds of problems. I try to isolate what it is that affects the problem and
> make the problem as simple as possible. Whatever is left must then at least
> be related to the cause.
>
> Best
>  W.
>
> --
> 
> Wolfgang Bangerth  email: bange...@colostate.edu
>
>www: http://www.math.colostate.edu/~bangerth/
>
> --
> The deal.II project is located at http://www.dealii.org/
> For mailing list/forum options, see https://groups.google.com/d/fo
> rum/dealii?hl=en
> --- You received this message because you are subscribed to a topic in the
> Google Groups "deal.II User Group" group.
> To unsubscribe from this topic, visit https://groups.google.com/d/to
> pic/dealii/5hC7jODg-7k/unsubscribe.
> To unsubscribe from this group and all its topics, send an email to
> dealii+unsubscr...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.
>

-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.


[deal.II] Solving time dependent heat equation with MPI (PETsc)

2017-10-10 Thread Mark Ma

Dear experienced deal.II users and developers, 

I want to solve a heat equation in the time domain with distributed memory 
using MPI, but the results are incorrect. In order to do so, I reference 
tutorial step-23 for time updating method and step-40 for implementing MPI. 
May I ask whether my boundary condition is right or not? Should we do 
compress() after apply_boundary_values()? Thanks in advance!

With best regards,
Mark





In PETsc or Trilinos, how could we set the boundary conditions and initial 
conditions?


Since old_solution_T would be read in other places (not shown here), it is 
initialized with ghost cells; while old_solutuon_T_cal is only used for 
written, it does not have ghost cells. see code as following,
//
old_solution_T_cal.reinit (locally_owned_dofs,mpi_communicator);
old_solution_T.reinit 
(locally_owned_dofs,locally_relevant_dofs,mpi_communicator);
//
..


  template 
  void ThermalDiffusion::run ()
  {


setup_system();
assemble_system();


VectorTools::project (dof_handler, constraints, QGauss(degree),
  InitialValues_T(),
  old_solution_T_cal);

old_solution_T = old_solution_T_cal;

LA::MPI::Vector tmp (locally_owned_dofs,mpi_communicator);
LA::MPI::Vector forcing_terms (locally_owned_dofs,mpi_communicator);

for (timestep_number=1, time=time_step;
 time<= global_simulation_end_time;
 time+=time_step, ++timestep_number)
  {
pcout << "Time step " << timestep_number
  << " at t=" << time
  << std::endl;


//---
//run to solve T
//---

  //
  //time dependent

//assign right hand side
mass_matrix_T.vmult (system_rhs, old_solution_T_cal);

laplace_matrix_T.vmult (tmp, old_solution_T_cal);
system_rhs.add (-time_step * (1-theta), tmp);

assemble_rhs_T (time);
forcing_terms = dynamic_rhs_T;
forcing_terms *= time_step * theta;

assemble_rhs_T (time - time_step); 
tmp = dynamic_rhs_T;
forcing_terms.add (time_step*(1-theta),tmp);
system_rhs.add (1,forcing_terms);

//assign system matrix
system_matrix.copy_from (mass_matrix_T);
system_matrix.add (time_step * theta, laplace_matrix_T);



   

















* {  BoundaryValues_Temperature 
boundary_values_function_T;  //boundary_values_function.set_time 
(time);  std::map 
boundary_values_T;  VectorTools::interpolate_boundary_values 
(dof_handler,
BOUNDARY_NUM,
boundary_values_function_T,
boundary_values_T);  MatrixTools::apply_boundary_values 
(boundary_values_T,  
system_matrix,  
solution_T,  
system_rhs,  false);}*
solve_T ();

if (Utilities::MPI::n_mpi_processes(mpi_communicator) <= 32)
{
TimerOutput::Scope t(computing_timer,"output");
output_results ();
}

computing_timer.print_summary ();
computing_timer.reset();

pcout << std::endl;


old_solution_T = solution_T;


  }
  }
}

//definition of solve_T() 

  template 
  void ThermalDiffusion::solve_T ()
  {
TimerOutput::Scope t(computing_timer,"solve_T");
LA::MPI::Vector completely_distributed_solution 
(locally_owned_dofs,mpi_communicator);
SolverControl   solver_control (dof_handler.n_dofs(),1e-12);

#ifdef USE_PETSC_LA
LA::SolverCG solver(solver_control,mpi_communicator);
#else
LA::SolverCG solver(solver_control);
#endif

LA::MPI::PreconditionAMG preconditioner;
LA::MPI::PreconditionAMG::AdditionalDatadata;

#ifdef USE_PETSC_LA
data.symmetric_operator = true;
#else
/* Trilinos defaults are good */
#endif

preconditioner.initialize(system_matrix,data);

solver.solve 
(system_matrix,completely_distributed_solution,system_rhs,preconditioner);



pcout << "  Solved in " << solver_control.last_step()
  << "  iterations."<< std::endl;

constraints.distribute (solution_T);

solution_T = completely_distributed_solution;

pcout << "success...solve_T()..." <