[aspect-devel] The output of nonadiabatic_temperature

Juliane Dannberg dannberg at gfz-potsdam.de
Fri Oct 14 10:34:49 PDT 2016


Hi Nan,

I admit that I haven't looked over all of your code in detail, but one 
point I saw is that you have to use this->get_geometry_model() and not 
this->get_geometry().
Hope that helps!

Best,
Juliane


On 10/14/2016 05:11 AM, Nan Zhang wrote:
> Hi all,
>
> In order to output the excess temperature, I changed 
> the nonadiabatic_temperature.cc like:
> ********************************************************************************************
> ...
> template <int dim>
>       void
>       NonadiabaticTemperature<dim>::
>       compute_derived_quantities_vector (const 
> std::vector<Vector<double> >  &uh,
>                                          const 
> std::vector<std::vector<Tensor<1,dim> > > &,
>                                          const 
> std::vector<std::vector<Tensor<2,dim> > > &,
>                                          const std::vector<Point<dim> > &,
>                                          const std::vector<Point<dim> 
> >  &evaluation_points,
>  std::vector<Vector<double> >  &computed_quantities) const
>       {
>         const unsigned int n_quadrature_points = uh.size();
>         Assert (computed_quantities.size() == n_quadrature_points,   
>  ExcInternalError());
>         Assert (computed_quantities[0].size() == 1,             
> ExcInternalError());
>         Assert (uh[0].size() == this->introspection().n_components, 
> ExcInternalError());
>
>         for (unsigned int q=0; q<n_quadrature_points; ++q)
>           {
>             const double temperature = 
> uh[q][this->introspection().component_indices.temperature];
>
> //computed_quantities[q](0) = temperature - 
> this->get_adiabatic_conditions().temperature(evaluation_points[q]);
>             double depth = 
> this->get_geometry().depth(evaluation_points[q]); //position???
>             const unsigned int idx = static_cast<unsigned 
> int>((ave_temp.size()-1) * depth / 
> this->get_geometry_model().maximal_depth());
>             double delta_temperature = temperature-ave_temp[idx];
> computed_quantities[q](0) = delta_temperature;
>           }
>       }
>
> ...
> **************************************************************************************
>
> A major problem is the data structure for the 
> this->get_geometry().depth(?). The cell type evaluation_points[q] does 
> not work. From the library definition, this->get_geometry().depth(?) 
> needs a position type array inside. How to give such an array to get 
> the depth information. I update the ave_temp array outside of 
> this subroutine with 50 layers.
>
> Bests,
> Nan
>
> On Wed, Oct 12, 2016 at 1:50 AM, Juliane Dannberg 
> <dannberg at gfz-potsdam.de <mailto:dannberg at gfz-potsdam.de>> wrote:
>
>     Hi Nan,
>
>     I haven't found a good way to plot the depth average files in
>     Paraview, and as far as I know we also don't have a visualization
>     postprocessor that does exactly what you want to do (if anyone has
>     a good idea please speak up!)
>     When I used the "depth average" postprocessor I always used pyplot
>     for visualization, so if you already have all of your output as
>     .csv you could write a python script to make your plot.
>
>     If you want to use Paraview, you would probably have to write your
>     own postprocessor in Aspect that outputs exactly what you want
>     (and we would be happy to see a pull request for that!).
>     Alternatively, you could write your own adiabatic conditions
>     plugin that computes the average temperature for each depth and
>     updates it every time step and uses that as adiabatic temperature
>     (and then you could use the "nonadiabatic temperature"
>     postprocessor). That second option is more of a hack, and would
>     only work if you don't need the adiabatic temperature for anything
>     else in your model, because in this case your adiabatic conditions
>     would contain something that's not really an adiabatic profile any
>     more, but it would probably be easier to implement.
>
>     Best,
>     Juliane
>
>
>     On 10/11/2016 11:05 AM, Nan Zhang wrote:
>>     Hi Juliane,
>>
>>     Thanks for your suggestion. What I really want is the excess
>>     temperature, the T subtracted from the average temperature on
>>     depth. I have done the post-processing with the .csv format,
>>     because Paraview converts the .vtu file to .csv file.
>>
>>     My problem is the csv file is not good for paraview plotting
>>     anymore. Paraview plots csv file with table to point. There is no
>>     any smooth/grid between points. The visualization is very bad.
>>     Have you guys done any conversion from csv back to vtu before??
>>
>>     Cheers,
>>     Nan
>>
>>     On Wed, Oct 12, 2016 at 12:53 AM, Juliane Dannberg
>>     <dannberg at gfz-potsdam.de <mailto:dannberg at gfz-potsdam.de>> wrote:
>>
>>         Hi Nan,
>>
>>         there are different options for outputting the temperature
>>         difference between plume and background mantle, and which one
>>         you want to use depends on what exactly you want to know:
>>
>>         If you use the "nonadiabatic temperature", you will get the
>>         excess temperature of the plume with respect to an adiabatic
>>         mantle temperature profile. If you are worried your adiabat
>>         will change over time, there is an update() function in the
>>         interface of the adiabatic conditions, so in principle you
>>         could use an existing model for the adiabatic conditions and
>>         implement this function, and then your adiabatic profile
>>         would be updated every time step.
>>
>>         Alternatively, if you want to compare the plume temperature
>>         to the current average mantle temperature at a given depth,
>>         you can use Aspect's "depth average" postprocessor. It
>>         computes depth averaged quantities (including the
>>         temperature) and writes them into a separate output file.
>>
>>         Best,
>>         Juliane
>>
>>
>>         On 10/11/2016 04:23 AM, Nan Zhang wrote:
>>>         Hi all,
>>>
>>>         I am doing a compressible model and try to plot the plume
>>>         structure. The plume is defined by the high temperature
>>>         difference from current average temperature at every depth.
>>>         I see ASPECT has an output "nonadiabatic_temperature". I
>>>         wonder if this output serves my purpose?
>>>
>>>         What I concerned is if the adiabatic temperature is not
>>>         exactly the same as my average temperature at every depth.
>>>         When I set up my model, I initialize the adiabatic
>>>         temperature profile with specific parameters. In theory, it
>>>         should be the same as the average temperature at every
>>>         depth. But, after billion year calculation, the average
>>>         temperature at every depth deviates away from the initial
>>>         adiabatic temperature profile.
>>>
>>>         So, I wonder if there is an output
>>>         "average_subtracted_temperature" in ASPECT? If so, it could
>>>         also serve the plume in the incompressible convection model.
>>>
>>>         Bests,
>>>         Nan
>>>
>>>
>>>         _______________________________________________
>>>         Aspect-devel mailing list
>>>         Aspect-devel at geodynamics.org
>>>         <mailto:Aspect-devel at geodynamics.org>
>>>         http://lists.geodynamics.org/cgi-bin/mailman/listinfo/aspect-devel
>>>         <http://lists.geodynamics.org/cgi-bin/mailman/listinfo/aspect-devel>
>>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.geodynamics.org/pipermail/aspect-devel/attachments/20161014/6c2bec97/attachment-0001.html>


More information about the Aspect-devel mailing list