[CIG-SHORT] Adding fields to HDF5/Xdmf files

[Brad Aagaard]

On 12/04/2013 02:10 PM, Ma, Xiao wrote:
> Hi Brad,
> This what I am write about the script , but it doesn't work,
> import numpy
> import h5py
>
> # Load in change in tractions from coseismic simulation
> h5 = h5py.File("output/stress.h5", 'r', driver="sec2")
>
> vertices = h5['geometry/vertices'][:]
> #tractions_change = h5['vertex_fields/traction_change'][0,:,:]
> stress=h5['cell_fields/stress'][:]
> #h5.close()
> pstress=(stress[:,:,0]+stress[:,:,1])/2+numpy.sqrt(((stress[:,:,0]-stress[:,:,1])/2)**2+(stress[:,:,2])**2)
> field = h5.create_dataset('vertex_fields/pstress',data=pstress[:,:])
                              ^^^^^^^^^^^^^ This should be cell_fields 
because it is a field over the cells.

> field.attrs['vector_field_type']='tensor'
                                     ^^^^^^^ You may need to change this 
to 'other' because you are in 2-D and not 3-D.

> h5.close()

You are missing the pylith_genxdmf part (that can be run from the 
command line).

NOTE 1: PLEASE ALWAYS SEND QUESTIONS TO THE EMAIL LIST 
cig-short at geodynamics.org.

NOTE 2: You should always include any error messages and detailed 
explanations of what is wrong. When you say "it doesn't work", I have to 
guess at what is happening.

Regards,
Brad



>
> best
> Xiao
> ________________________________________
> From: Brad Aagaard [baagaard at usgs.gov]
> Sent: Wednesday, December 04, 2013 15:47
> To: Ma, Xiao; cig-short at geodynamics.org
> Subject: Adding fields to HDF5/Xdmf files
>
> Xiao,
>
> You can add new fields to the HDF5 file using the create_dataset
> function in the h5py module. To update the Xdmf file, use the
> pylith_genxdmf script supplied with PyLith:
>
> For example, in a Python script I do something like:
>
> field = h5.create_dataset('vertex_fields/pstress', data=pstress[:,:,:])
> field.attrs['vector_field_type'] = 'tensor'
> h5.close()
>
> # Regenerate XDMF file.
> import os
> os.system('pylith_genxdmf --file=statevars.h5')
>
>
> Regards,
> Brad
>
>
> On 12/04/2013 01:26 PM, Ma, Xiao wrote:
>> Hi Brad,
>> What I mean is that I want to use the principle values to generate the Xdmf file , to visualize in the paraview.
>> Best,
>> Xiao
>> ________________________________________
>> From: Brad Aagaard [baagaard at usgs.gov]
>> Sent: Wednesday, December 04, 2013 14:57
>> To: cig-short at geodynamics.org; Ma, Xiao
>> Subject: Re: [CIG-SHORT] (no subject)
>>
>> Xiao,
>>
>> The stresses are computed at each of the quadrature points within a
>> cell. This is true for all cell related data. By default the output
>> includes the values at all of the quadrature points. You can use the
>> CellFilterAvgMesh filter (see the tutorials in the manual and the
>> examples) to have PyLith average the values over the quadrature points
>> so you end up with the values (e.g., stress tensor) at just point per cell.
>>
>> Regards,
>> Brad
>>
>>
>> On 12/04/2013 12:45 PM, Ma, Xiao wrote:
>>> Hi Brad,
>>> I am able to use the h5py , put the stress into a variable called stress, however the shape of this stress variable is [500,10000,12], I know the first tow 500 is the total times steps, and 10000 is the total element numbers, however why there is 12 components for the stress  tensor?
>>> Best,
>>> Xiao
>>>
>>>
>>>
>>> _______________________________________________
>>> CIG-SHORT mailing list
>>> CIG-SHORT at geodynamics.org
>>> http://geodynamics.org/cgi-bin/mailman/listinfo/cig-short
>>>
>>
>