[aspect-devel] The output of nonadiabatic_temperature
Juliane Dannberg
dannberg at gfz-potsdam.de
Fri Oct 14 22:45:00 PDT 2016
Hi Nan,
can you send the error message you get when you try to compile your code?
I don't see why this->get_geometry_model().depth (evaluation_points[q])
shouldn't work, that is the same way we do it for example in the depth
visualization postprocessor (source/postprocess/visualization/depth.cc).
evaluation_points is a vector of points, so if you access one element of
the vector using evaluation_points[q] you get the type Point<dim>, which
is exactly the input you need for the depth(Point<dim>) function.
Best,
Juliane
Am 10/14/2016 um 9:53 PM schrieb Nan Zhang:
> Hi Juliane,
>
> Yes, you are correct. I should use
> this->get_geometry_model().depth(evaluation_points[q]);
>
> But, my problem is the same. The compiling report that the
> get_geometry_model().depth() uses Point<dim> not vector<Point<dim> >.
> My understanding is the evaluation_points[q] is for the quadrature
> point. While the this->get_geometry_model().depth() should use kind
> of position data type, which I don't know!?
>
> Nan
>
> On Sat, Oct 15, 2016 at 1:34 AM, Juliane Dannberg
> <dannberg at gfz-potsdam.de <mailto:dannberg at gfz-potsdam.de>> wrote:
>
> 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/b940f6c6/attachment-0001.html>
More information about the Aspect-devel
mailing list