Hi there,
I met a problem when trying to apply the coseismic slip using linear interpolation, as shown below, in Pylith 5 and also Pylith 4.2 for 3D viscoelastic case, using displacement.basis_order = 2, lagrange_multiplier_fault.basis_order = 2 and defaults.quadrature_order = 3.
pylithapp.problem.interfaces.fault_slabtop.eq_ruptures.rupture
origin_time = 0.0*year
db_auxiliary_field = spatialdata.spatialdb.SimpleDB
db_auxiliary_field.description = Coseismic slip
db_auxiliary_field.iohandler.filename = Element_reverse_linear_330.spatialdb
db_auxiliary_field.query_type = linear
When db_auxiliary_field.query_type = nearest, the coseismic simulation runs without any problem but when I try to use db_auxiliary_field.query_type = linear, error always occurs, as shown below.
[0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------
[0]PETSC ERROR: Error in external library
[0]PETSC ERROR: Could not find values for initiation_time at ( -247691 588111 -46561.8) in spatial database ‘Coseismic slip’.
I learned during workshop that the coseismic slip is better to be linearly interpolated on fault nodes to be more natural. The Element_reverse_linear_330.spatialdb file contains the location of element center of all triangles elements on the target fault interface and only within a small patch (a triangle area contacting 16 elements, for Green function purpose) final_slip_reverse is set as 1 and all other area has final_slip_reverse=0. I wonder how to successfully use db_auxiliary_field.query_type = linear, and is there anything I missed? Thank you in advance for your attention.
Best,
Qingjun Meng
In PyLith v3 and later, the default basis order is 1 for all fields. So the default basis order for the auxiliary field for the slip time history parameters (initiation time and final slip values for the default step function) will be 1. This means the field values will be at the vertices, not the cell centroids.
In most cases, we specify the spatial database on a grid of points that includes the rupture patch but may be at different points than the mesh. For example, we might use a uniform grid of points on the fault surface even though the mesh on the fault surface is triangles. As long as both are smooth and have sufficient resolution, the discretization error will be small.
Thank you for your response, Brad. So you suggest that I need to use SimpleGridDB format and project input slip points on evenly spaced grids in 3D, instead of using SimpleDB format with slip on element centers to do the linear interpolation? For a 2D problem, the Simple DB format seems work OK, for example the subduction-2d case. Maybe the 3D problem is more complex?
Qingjun
The philosophy behind a spatial database is to be able to specify fields (slip, material properties, etc) in the coordinate system and geometry of the data. We have spatial databases for analytical functions, data on a uniform grid, and data at arbitrary locations. Both SimpleGridDB and SimpleDB spatial databases can be used for uniformly gridded data. The SimpleGridDB is more efficient for queries, but is more difficult to setup for specifying slip on a complex fault surface in 3D.
I’m not sure if I fully understand your point. Previously I thought uniform grid should be equal dx, dy and dz, which means SimpleGridDB. But maybe you suggest a uniform grid of dx and dy (for 3D reverse fault plane), keeping z=0 with data-dim = 2, then the linear interpolation is done over a simple horizontal plane not necessary to consider about the depth info of fault interface, as long as the grid points are dense enough to keep the resolution.
Thanks,
Qingjun Meng
A SimpleGridDB is a logical grid (points do not have to be evenly spaced) but the grid must have locations that can be described by arrays of x, y, and z when working in 3D. A SimpleDB, on the other hand, can have points on a logical grid but the locations are only given by the coordinates of the points, making it easier to specify points on an arbitrarily oriented surface (perhaps nonplanar) in 3D.
Thank you, Brad. Yes, SimpleGridDB does not need to have uniform grid. But for SimpleDB, the uniform grid is only applied for x, y coordinates and z coordinate is the real depth (based on the fault interface 3D geometry) or could be to zero? I found for the 2D problem, coseismic slip .spatialdb can only have y coordinate and x coordinate is all zero. The same way could be applied to 3D problem, after I set data-dim=2 ,right? The linear interpolation is done after slip-impulse data points got projected onto a simple surface, which could be the z=0 surface.
Qingjun
For a SimpleDB, you can use a uniform grid for any coordinate. The implementation used makes no assumptions about the layout of the points other than the spatial dimension and the data dimension, which are specified in the header.
You can specify slip on a vertical or horizontal plane or on a plane that roughly matches the fault. The points in the mesh are projected onto the plane (when data-dim == 2) of the points in the spatial database and then the field is linear interpolated to the point of interest.
Thank you, Brad. That is much clearer to me now. I could select a horizontal plane for the subduction interface. I really appreciate your instruction.