Segments with 0 `burial_depth` still have elastic Okada deformation even when replaced by a mesh.

[brendanjmeade]

This feels strange. When a group of segments has a burial_depth of 0 there still seems to be elastic deformation around them associated with Okada dislocations even when the entire group of segments should be replaced by a .msh file. This is unintended behavior.

Notes:

• This is observed with burial_depths equal to zero. Adjacent burial_depths of 15 do not show this behavior. I don't know and have not experimented with other burial_depths.

Expected behavior:

• burial_depth should have no influence in this case because the elastic deformation from the rectangular Okada segments should be ignored when they are replaced with a TDE mesh.

Potential solutions:

• If statement to truly set partials for these segments to zero?
• Set locking_depth and burial_depth both equal to zero for segments replaced by meshes?
• Simply eliminate any reference to burial_depth because we don't have a modern use for it?

[jploveless]

I'm not sure that it's just a burial_depth issue. I'm checking the sum of partials for replaced segments, and only a handful are actually zero. In the Japan example, locking depths are indeed getting set to zero for replaced segments (np.max(segment.locking_depth[toggle_off]) = 0), where toggle_off is the composite conditional statement from zero_mesh_segment_locking_depth, and burial_depth is all zero for these segments as well, but somehow there's a non-zero elastic contribution being calculated.

This is my code for checking these:

# Get indices of toggled-off segments, copied from zero_mesh_locking_depth
toggle_off = np.where((segment.patch_flag != 0) & (segment.patch_file_name != 0) & (segment.patch_file_name <= len(meshes)))[0]
print(np.max(segment.locking_depth[toggle_off])) # Confirm zero LD
print(np.max(segment.burial_depth[toggle_off])) # Confirm zero BD
sum_partials = np.sum(operators.slip_rate_to_okada_to_velocities[:, slice(0, 3*segment.shape[0], 3)], axis=0) # Sum rows of strike-slip partials
print(np.where(sum_partials == 0)) # Segment indices of zero partial sums
print(toggle_off) # Segment indices of toggled-off segments
print(np.intersect1d(np.where(sum_partials == 0), toggle_off)) # Intersection of indices 

I'm not sure what I'm missing or getting wrong here, but I'll keep checking.

[brendanjmeade]

@jploveless I think I'm understanding this better now. In zero_mesh_segment_locking_depth(segment, meshes) we have the condition:

    toggle_off = np.where(
        (segment.patch_flag != 0)
        & (segment.patch_file_name != 0)
        & (segment.patch_file_name <= len(meshes))
    )[0]

For the case that I'm working on (Anatolia) toggle_off isn't toggling off the segments I expect and the reason for that is that there is only one patch_file_name and its name is 0. I think that's my big problem. Are you using 1 as the name of the first patch_file_name?

[jploveless]

Yes — I am using 1 as the minimum patch_file_name. I did this as a holdover from Matlab, but I think that we might be able to use patch_flag only. Using the file name was important when the edge constraints and smoothing weights were specified just as a vector in the .command file, but now that we use the mesh parameter file for those conditions, we could probably just rely on patch_flag to mark the replacements.

I just did a search through the Blocks repo, and it doesn't seem that Segment.patchFile is used for anything terribly important. I think I liked to use it with SegmentManager just for a visual representation of how segments were associated with meshes. I think this means that we could define:

toggle_off = np.where(segment.patch_flag != 0)[0]

and this would eliminate the zero-based vs. one-based indexing issues that seem to be causing the problem.

Getting back to your original suggestion, I do think it's important to set burial_depth to zero for replaced segments, or else the segment will have a non-zero down-dip width and therefore will have an elastic contribution.

[jploveless]

Murphy's law: I just received an email about Blocks that reminded me of use of the patch_file_name field: the Blocks function TogglePatches.m looks through the segment structure, and if there is a mesh that is listed in the .command file (or .mshp file), but none of its corresponding segments have their toggles turned on, this mesh will be ignored.

A use case for this is a multiple-mesh model where the .segment file is modified to turn off the segment flags for a given mesh, but the same list of meshes is still read in as specified in the .command or .mshp file.

[brendanjmeade]

Yes — I am using 1 as the minimum patch_file_name. I did this as a holdover from Matlab, but I think that we might be able to use patch_flag only. Using the file name was important when the edge constraints and smoothing weights were specified just as a vector in the .command file, but now that we use the mesh parameter file for those conditions, we could probably just rely on patch_flag to mark the replacements.

I like the separation between patch_flag and patch_file_name. They indicated different things, and so I think it's reasonable to treat them as separate variables. Future us will like this logical clarity!

I think that the "first" patch_file_name should be 0 because it would enable us to use patch_file_name as an index into the meshes list. I don't think we do that now, but that could be super useful in the future.

I just did a search through the Blocks repo, and it doesn't seem that Segment.patchFile is used for anything terribly important. I think I liked to use it with SegmentManager just for a visual representation of how segments were associated with meshes. I think this means that we could define:

toggle_off = np.where(segment.patch_flag != 0)[0]

I like this a lot and think it's the correct design decision. Feel free to push!

and this would eliminate the zero-based vs. one-based indexing issues that seem to be causing the problem.

This makes sense looking at the code. I think the only place we might have issues is in snap_segments where we see (segment.patch_file_name == i + 1) as a part of a conditional. I think we should agree on a clear numbering routine for mesh files. I think starting with 0 would be very pythonic. Happy to discuss this.

Getting back to your original suggestion, I do think it's important to set burial_depth to zero for replaced segments, or else the segment will have a non-zero down-dip width and therefore will have an elastic contribution.

This makes sense. An alternative could be to not touch the locking depths or burial depths at all. Instead, we could just have an if check for patch_flag at some layer into the Okada elastic calculation that would simply return zeros for the cases with where there is a mesh replacing a segment. One plus side of this would be that we would never be modifying the locking depth or burial depth information through a run...actually, I guess we still can with locking depth flags. One downside of this would be that it might require patch_flag at an awkward point on the Okada calculation stack that we might not want for a lot of other reasons. Still writing this down but I think zeroing out the locking depths and burial depths for these segments might be the more reasonable thing to do because it would produce an output segment file that's immediately useful for forward calculations.

Happy to consider other options!

[jploveless]

I like the separation between patch_flag and patch_file_name. They indicated different things, and so I think it's reasonable to treat them as separate variables. Future us will like this logical clarity!

I think that the "first" patch_file_name should be 0 because it would enable us to use patch_file_name as an index into the meshes list. I don't think we do that now, but that could be super useful in the future.

I think this would be fine, but I'm a little concerned about the segments that don't correspond to any patch being assigned a patch_file_name of 0 by default. Should we set those to something like NaN or -999 or just -1? I think this change should be implemented in celeri_ui as well as default behavior when generating a new segment.

I like this a lot and think it's the correct design decision. Feel free to push!

I'll hold off on the toggle testing change until we have a better sense of how to index patch_file_name. If we decide on a negative number for all segments that don't have a mesh, we could keep the composite conditional, just changing it to find segment.patch_file_name >= 0 instead.

This makes sense looking at the code. I think the only place we might have issues is in snap_segments where we see (segment.patch_file_name == i + 1) as a part of a conditional. I think we should agree on a clear numbering routine for mesh files. I think starting with 0 would be very pythonic. Happy to discuss this.

Agreed. This is one place where it's important to keep the reference to the mesh file name. It would be straightforward to just edit snap_segments to correspond to zero-based indexing.

This makes sense. An alternative could be to not touch the locking depths or burial depths at all. Instead, we could just have an if check for patch_flag at some layer into the Okada elastic calculation that would simply return zeros for the cases with where there is a mesh replacing a segment. One plus side of this would be that we would never be modifying the locking depth or burial depth information through a run...actually, I guess we still can with locking depth flags. One downside of this would be that it might require patch_flag at an awkward point on the Okada calculation stack that we might not want for a lot of other reasons. Still writing this down but I think zeroing out the locking depths and burial depths for these segments might be the more reasonable thing to do because it would produce an output segment file that's immediately useful for forward calculations.

I do like the idea of formally zeroing out the locking and burial depths for replaced segments, so that the saved versions of the segment file reflect those changes. I think that any changes we make to the geometry within a celeri run should be preserved as part of the output. As we move toward using more and more meshes, it may make sense to only run the Okada routine for those segments with non-zero locking and burial depths.

[brendanjmeade]

@jploveless Thanks for all of these ideas! Agreed on all points. Below I've tried to summarize and convert to a finite work list. Let me know what you think

Things that can happen now

• [x] Set locking_depth and burial_depth to 0 for all segments with patch_flag == 1 in zero_mesh_segment_locking_depth()
• [x] Adjust indexing in snap_segments to work with 0-based indexing
• [x] In celeri_ui set default patch_file_name to -1.
• [x] In zero_mesh_segment_locking_depth modify conditional to be:

  toggle_off = np.where(
  (segment.patch_flag == 1) # Is this sufficient?
  & (segment.patch_file_name >= 0) # new condition
  & (segment.patch_file_name <= len(meshes))
  )[0]

Things to consider

• For existing *_segment.csv files: write a script to update to new standard:
  • Read file: *_segment.csv
  • Find all patch_flag == 1 segments
    • Find min and values of patch_file_name for each of these segments
      • If min patch_file_name == 0 do not update
      • If min patch_file_name == 1 decrement all patch_file_name values by 1 so that they now start at 0
  • For all segments with patch_flag != 1 set patch_file_name to -1
  • Write file: *_segment.csv
• Run this on all *_segment.csv files in celeri/data/segment

[jploveless]

https://github.com/brendanjmeade/celeri/commit/f5fbf74b8ff7bb4a8395913652b620f4a23c216e addresses the minor code changes (items 1, 2, 4 of the checklist).

I like the workflow for modifying all existing segment files but I haven't done that yet.

[brendanjmeade]

Wow @jploveless! I'll do the celeri_ui bit (fingers crossed!) and script to update existing segment files.

[brendanjmeade]

celeri_ui default values for patch_file_name set to -1 in https://github.com/brendanjmeade/celeri_ui/commit/cd9f392639b7040b0c9ba9d452e4fe1a4df3d6a3

[brendanjmeade]

Completed with https://github.com/brendanjmeade/celeri/commit/745da012b032159afae6f311ed2a2dc377addccf

• Script for updating *_segment.csv files /utilities/celeri_patch_file_name_adjust.ipynb
• Updated *_segment.csv files in /data/segment