Smoothing two-stack microstructure - jutting nodes after smoothing

[sgbaird]

Hi @siddharth-maddali,

Thanks for your help so far. I've been plotting the grains after doing h-smoothing and have been getting some strange results. See for example:

I also get:

109030 of 800612 nodes not smoothed.

which may be related to the quad points being too close together as you mentioned before, or it could be an issue of having a thickness of only 2 voxels in the z-dimension (I have EBSD data from both sides of a sample with a large proportion of through-thickness grains).

Using the sample data produces expected results:

If I apply DREAM.3D's laplacian smoothing filter before doing the hierarchical smoothing, I get the following, where on further inspection, it seems the triangles that "jut" leave concave "divots" in the grain.

with the following workspace:

Name                   Size                   Bytes  Class                        Attributes

  GrainToPlot            1x1                        8  double                                 
  fl               1438608x2                 23017728  double                                 
  ngrains             3348x1                    26784  double                                 
  s1                     1x1                        8  matlab.graphics.axis.Axes              
  s2                     1x1                        8  matlab.graphics.axis.Axes              
  tri              1438608x3                 34526592  double                                 
  trisub               852x3                    20448  double                                 
  xdat                   3x800607            19214568  double                                 
  xsmooth                3x800607            19214568  double     

One round-about way I might deal with this is turning the 2-stack into a 4-stack by copying the stacks above and below as appropriate and cropping it afterwards, but I'm worried this would have the artificial effect of the grain boundary always being normal to the surface when it's close to the surface.

Another similar approach would be breaking up the interior into multiple voxels, though I'm not sure how I would go about this.

Do you think this is an inherent limitation of the method for this unique use case?

[siddharth-maddali]

@sgbaird it's difficult to tell without more information on your volume. How are the nodes labeled in your example? What are the ntype labels of the nodes that "shoot out into space" after you run the smoothing?

This might have something to do with what you're seeing. h-smooth believes whatever mesh you give it in the tri variable without checking to see if the topology is messed up. I suspect some of the "jutting" nodes are labeled as interior points (2 or 12), which have the greatest freedom to move around in space.

Since you're dealing with a 3-layer (2-voxel) structure, you should try this: artificially re-label the nodes on the "rim" in all 3 layers as a triple line (3, 13), possibly bound by quad points (4, 14) if known in advance, and then run h-smooth. By this I mean use a new, custom ntype array. Note that the interior nodes in each layer can be labeled 2/12 as usual, as long as they're surrounded in the plane by nodes of type 3/13 or 4/14. If you know exactly which points along your rim are of type 4 (quad points), you should add that in artificially. Doing this hack will effectively tell h-smooth that your volume has 2 grains stacked one on top of the other, each of single-voxel thickness.

Back in the day I didn't bother adding a 2D smooth capability to h-smooth, because you can easily "hack" the 3D smooth for 2D meshes with the correct re-labeling of the nodes. Proper smoothing really depends on the tri-mesh you give it (remember, this remains unchanged before and after smoothing). If the mesh was defined poorly, the smoothing will proceed unpredictably. This problem goes back to DEREAM.3D's QuickMesh, which behaves funny sometimes.

Let me know if this hack works...

[sgbaird]

Hi @siddharth-maddali,

Thanks for your suggestions!
I modified the Test function to plot the node types (see attached)
run.txt
Test.txt

Isoview:

Sideview:

Another example:

Based on this it seems pretty apparent that jutting nodes are interior nodes (blue) as you suspected.
Partially for my own reference and partially for others, I'm copying the node definitions here:

Id Value	Node Type
2	Normal Vertex
3	Triple Line
4	Quadruple Point
12	Normal Vertex on the outer surface
13	Triple Line on the outer surface
14	Quadruple Point on the outer surface

Also, since the first grain was a pretty good candidate to show what I'm trying to accomplish on a higher level, here is the side-view "interpolation" that I'm essentially hoping to get in the end by whatever means necessary:

One thing I didn't explain (or I may not have understood your comment about 2 grains being stacked on top of each other) is that the top and bottom layers should generally be the same grain (i.e. a columnar or "through-thickness" grain).
If h-smooth interprets the volume as "2 grains stacked one on top of the other, each of single-voxel thickness", would it still tend towards the interpolation I'm hoping to accomplish above?

[sgbaird]

Also, to give a little more context, DREAM.3D's Laplacian smoothing filter starts to do what I'm hoping for, but despite trying different parameters the outer surfaces either always retained their "blocky" appearance or made strange distortions to the outer surfaces (latter not shown here). This (and a post on the DREAM.3D google groups forum) were what originally led me to look at using h-smooth.
Isoview (hsmoothed version on right):

Topview:

[siddharth-maddali]

Hi @sgbaird ,
I think I understand your data a little better now. Based on the side view you've shown, I take back what I said earlier about pretending to have two grains, one stacked on the other. That was the wrong visual to have it seems, sorry for that!

There might be a problem because of the size of the Z-axis separation of the layers, relative to the mesh size in the XY plane (assuming the isosurfaces you've shown have aspect ratio 1:1:1). You might try artifically compressing the Z-dimension of the grain to be comparable to the size of the XY discretization (i.e., the size of the mesh element). This will involve changing the Z components of your nodes. Then attempt the smoothing, and finally apply the inverse Z-scaling operation on the smooth nodes. By doing this, you might stop the nodes jutting out like that. I can't guarantee that you'll get the smooth variation from layer to layer that you're looking for, but you might reduce the "steppiness".

If the nodes still jut out, then I'm out of ideas! Other than blaming it on possibly poor meshing, because I know for sure that if the mesh is topologically sound, h-smooth will work! I am, however, curious. Can you share a data file for that particular grain (nodes, tri, ntype, fl) that I can load in Matlab? I'm quite busy these days, but I'll play with it in my free time and see what I can glean from it.

[imikejackson]

Just to comment that we started to integrate the H-smooth into DREAM3D here but ran into problems that stem from a commit that we put into the QuickMesh to try to "fix" some bad nodes (Commit 1ed41586a0bf7e9c528a11dd65c2530aa1d7682d) which I think we got as far as figuring out that this caused some of the issues with the H-Smooth algorithm. Unfortunately I don't have any kind of timeline as to when the code will be completed.

[sgbaird]

Apparently @jercgreen6 ran into some problems with nodes jutting out, even though the microstructure he was using is fairly straightforward (cube with "normal" looking grains). Additionally, he mentioned that he's only seeing the issue with the C++ code, but not the native (slow) MATLAB version. @jercgreen6, could you comment more on this?

[jercgreen6]

Initially I set up the MATLAB version of Hierarchical Smooth to verify that the code would do what I wanted it to do. Other than computation time the results were ideal so I set up the C++ version. I modified my MATLAB script which I use to run Hierarchical Smooth to match the input arguments and remove the buffer elements. However, I get very different results between the two codes for the same microstructure. The C++ version produced spikes in the microstructure similar to what was shown by @sgbaird.

I also get the same warning:
WARNING xxxx of xxxx nodes not smoothed.

This has occurred for all of the microstructures that I have tried. Microstructures were produced using Dream3D.

Figure 1: Smoothed microstructures using the C++ version of the code

Figure 2: Smoothed microstructure using the MATLAB version of the code

Included in the zip folder is the data that I used as the input to Hierarchical Smooth.
M.VerticesSmooth = HierarchicalSmoothMatlab( tri, M.Vertices.', fl, nodetype).';
M15Data.zip

@siddharth-maddali Do you have any ideas on what is going wrong here? Do you have suggestions on how to address this issue?

Thanks!

[sgbaird]

@jercgreen6, nice explanation. I wonder if you get the same kind of results using ex1 and ex2 from the examples in the repository using these instructions. If it runs OK for the example microstructures with the C++ code, it makes me think maybe it's an issue with the export/conversion of labels .. Might be worth checking to see if there are any blaring differences between the text files in the example repositoyr and M15Data. Also, did you remove the buffer mesh elements on the faces of the box:

f = find( any( fl > 0, 2 ) );
tri = tri( f, : );
fl = fl( f, : );

before running it through HierarchicalSmoothMatlab()? This might not affect anything (the documentation suggests it only affects the runtime), but thought I'd check.

All-in-all a strange issue and the inconsistency between versions suggests a possible bug. How urgent is it (days, weeks, etc.)?

Also, @siddharth-maddali, sorry for not including a dataset before. Hopefully the one @jercgreen6 included will be enough to resolve the issue. Thanks in advance for the help! Let us know when you get this and what your thoughts are

[jercgreen6]

For the examples I still get the error about some of the nodes not smoothing in the C++ code (they work fine for the MATLAB code). This would suggest that it is either how I am interfacing with the code or an issue with the code itself. I am in the troubleshooting phase to try an figure out what is going on.

I do remove the buffer elements before running the code. I tried running it without removing the buffer elements and more of the nodes smoothed but there were still a large amount which did not.

A solution within 2 weeks would be ideal but not required.

[sgbaird]

@jercgreen6 do you still get the "jutting nodes" with the examples via C++?

@siddharth-maddali had mentioned to me in #9 (comment) about why some number of nodes aren't smoothed, but this didn't seem to be a major problem.

[sgbaird]

It would appear there's also a difference in indexing the triangulation between the MATLAB and C++ versions.

For the MATLAB version, based on the documentation, I think it's:
​tri = ​1​ ​+​ load( ​'​../../examples/ex2/SharedTriList.txt​'​ );

Whereas for C++ MEX version I'm guessing it's just:
​tri = load( ​'​../../examples/ex2/SharedTriList.txt​'​ );

Just trying to think through any possibilities other than an error in the C++ code which will probably be more painful to figure out.

[jercgreen6]

@sgbaird
So far I have been only able to observe the effect of the C++ code on the example files. The code which I use in conjunction with Hierarchical Smooth doesn't work well with the MATLAB version of Hierarchical Smooth. As far as I can tell this is due to several changes in Dream3D output files from when HierarchicalSmooth was created.

• Dream3D no longer creates microstructures containing buffer elements. These are not considered when using Hierarchical Smooth but buffer elements need to be removed from the dataset to process correctly with my code.
• Dream3D doesn't distinguish a grain neighboring empty space between the sides and top (only uses -1, no 0s). 0's need to be replaced by -1's.
• With the MATLAB code I am encountering an unknown issue due to the existence of several grains having grainID 0 (breaks MATLAB's indexing rules). I tried adding 1 to the entire array to just to see what would happen but it still wouldn't go through.

As far as your comments on indexing being the issue. I tried running the code with and without adding 1 to the tri variable. The results are shown here:

Figure 1: Example 1 where 1 was not added to the tri variable. Here nodes seem to be pulled in a single direction.

Alternatively, the closest that I have been able to get to things actually working are shown here. Previously I was using the fscanf command instead of the load command. I don't know why this would make a strong difference bet it seemed to help.

Figure 2: Example1 where load was used instead of fscanf to initialize variables. Still contains the ripple effect.

At this point all that I know is that the smoothed vertices that I am getting out of Hierarchical Smooth do not match what is recorded as the example output. I'm not sure if that is to be expected using the Hierarchical Smooth algorithm. About 25% of vertex locations are the same (determined using nnz(MySmoothedVertices == ExampleSmoothVertices)) but 75% are different. @siddharth-maddali After running Hierarchical Smooth should my output files exactly match the example smooth file?

[sgbaird]

@jercgreen6, looks like you're getting closer. Assuming you're using the latest version now, perhaps try using a 2018 version of DREAM.3D? From what I'm seeing, it doesn't seem like an internal error in the C++ code so much as a version incompatibility or maybe a data processing issue. I think a MWE for the DREAM.3D output process and HierarchicalSmooth which reproduces some of the plots you've been showing will make it easier for me to do some hands-on troubleshooting (if you have time to put a MWE together). I also think it's worth asking, have you tried any of the built-in DREAM.3D smoothing filters? (esp. Laplacian).

[siddharth-maddali]

@jercgreen6 @sgbaird Thanks for letting me know. I'll get back to you as soon as I can.

[siddharth-maddali]

@jercgreen6 ,

1. I don't think the Matlab code warns you that it's ignoring nodes, even though it does in the background. The printed message is a feature of the C++ code only, to warn you that your mesh might not be resolving your microstructure finely enough in some regions of your sample. I ran the C++ code on both example data sets and only ~30 out of 159433 nodes were not smoothed. This is reasonable because the data set has some very small grains. BTW, I used the C++ code wrapped in an Octave function, but the Matlab wrapper should act no differently.

2. It appears that your input M.Vertices has already been smoothed somewhat. Are you running HierarchicalSmooth on nodes that have already been already HierarchicalSmoothed? Or smoothed by some other process? If so, it might very well misbehave the second time round. HierarchicalSmooth makes use of the rectangular pixelation to constrain the smooth surface. It finds the best minimum-energy smooth surface that can 'slide' into the thin volume defined by a sheet (or line) of pixels. If this pixelation is gone, as it will be after one round of smoothing, then HierarchicalSmooth behaves in strange ways, and that could explain the jutting nodes.

3. I usually don't bother smoothing the flat cross-sections of the grains at the faces of the box, because these are not real grain boundaries anyway. As @sgbaird mentioned above, it is usually best to ignore these nodes while smoothing, and include them later in your plot.
   

I hope this helps! I need to look into why HierarchicalSmooth sometimes misbehaves for exactly flat surfaces, as it seems to be doing.

[siddharth-maddali]

@jercgreen6 I haven't worked with Dream.3d in a very long time, so I don't know about new changes to the indexing or labeling conventions. The code snippet above will give you the smoothed microstructure that you see in the animation at the top of README.md.

[siddharth-maddali]

@jercgreen6 Here is Octave code that deals with example 1, in which I've called the C++ routines in a wrapper. I'm not seeing significantly jutting nodes.

tri = 1 + load( '../examples/ex1/SharedTriList.txt' );
P = load( '../examples/ex1/SharedVertexList.txt' )';
fl = load( '../examples/ex1/FaceLabels.txt' );
nt = load( '../examples/ex1/NodeType.txt' );

% first, removing top/bottom buffers and keeping side walls
f = find( all( fl >= 0, 2 ) );
fl_side = fl( f, : );
tri_side = tri( f, : );


% Doing HierarchicalSmooth
tic; 
Ps = HierarchicalSmoothOctave( tri_side, P, fl_side, nt );
t = toc; 
printf( 'Runtime = %e seconds. \n', t );

% Creating a random color selection for grains in surface plot
col = max( fl_side, [], 2 );
col_u = unique( col );
col_u = col_u( randperm( numel( col_u ) ) );
[ ~, col ] = ismember( col, col_u );
save( 'example1_demo.h5', '-hdf5' );

[sgbaird]

I finally got around to trying out the C++/MATLAB version again using the ex2 data, and got the following:

While inside the Src/Cpp/ folder (same directory for everything else)

folder = fullfile('..', '..', 'examples', 'ex2');

xdat = load( fullfile(folder, 'SharedVertexList.txt') ).';

tri = 1 + load( fullfile(folder, 'SharedTriList.txt') );

fl = load( fullfile(folder, 'FaceLabels.txt') );

f = find( any( fl==-1 | fl==0, 2 ) );
fl( f, : ) = [];
tri( f, : ) = [];

ntype = load( fullfile(folder, 'NodeType.txt') );

whos

xsmooth = HierarchicalSmoothMatlab( tri, xdat, fl, ntype );

t = num2cell(xsmooth, 2);
trisurf(tri, t{:}, 'EdgeColor', 'none')
axis equal

I compiled it on Linux with a modified Makefile (just the CC= ... line changed)

CC= /apps/gcc/6.3.0/bin/g++

in order to compile with gcc 6.3. I modified CreateMatlabMex.m as follows and used MATLAB R2019b (on my compute cluster, module load matlab/r2019b):

mex -glnxa64 ...
    -I./Eigen ...
    -I./libigl/include ...
    -I/usr/include ...
    -I. ...
    -v GCC='/apps/gcc/6.3.0/bin/g++' ...
    HierarchicalSmoothMatlab.cpp libhsmooth.a

I had the most recent versions of Eigen and libigl cloned into my folder via:

git clone https://gitlab.com/libeigen/eigen.git

(per the instructions at Eigen)
and

git clone https://github.com/libigl/libigl.git

At first glance, the plot looks OK to me. Nothing really jumps out as being unphysical. @jercgreen6, can you run the MATLAB code (first code snippet) and post the plot that you get? Curious if you're getting the same output. If you have some code you want me to try, I can give that a go.

[jercgreen6]

@sgbaird
I get the same plot as you do for the example 2 data.

However, if I use the exact same code for files that I made with Dream3D I still encounter issues. For half of the microstructures I get some type of error which stops the code. The code runs to completion for the other microstructures but they come out spiky. I think the changes may be caused by variations in the Dream3D output files. The hierarchical smooth example files don't work with most of the code that I use for other Dream3D data and the Dream3D data that I make doesn't seem to be working well in Hierarchical Smooth.

I uploaded a couple .txt files that I used here:
https://drive.google.com/drive/folders/1k2RG4dHYZ9n9hRBVdBHsWJNHaoRQqQw9?usp=sharing
(One segmentation error, and the spiky microstructure shown above)

@sgbarid Do you get the same output that I do for the spiky microstructure?

[sgbaird]

@jercgreen6 try swapping the columns in FaceLabels.txt in the "spiky" dataset immediately after import using fl = fliplr(fl);

[sgbaird]

AKA

folder = fullfile('..', '..', 'examples', 'ex2');

xdat = load( fullfile(folder, 'SharedVertexList.txt') ).';

tri = 1 + load( fullfile(folder, 'SharedTriList.txt') );

fl = load( fullfile(folder, 'FaceLabels.txt') );
fl = fliplr(fl);

f = find( any( fl==-1 | fl==0, 2 ) );
fl( f, : ) = [];
tri( f, : ) = [];

ntype = load( fullfile(folder, 'NodeType.txt') );

whos

xsmooth = HierarchicalSmoothMatlab( tri, xdat, fl, ntype );

t = num2cell(xsmooth, 2);
trisurf(tri, t{:}, 'EdgeColor', 'none')
axis equal

FaceLabels.txt:

ex2	"spiky"
421, -1	-1, 23
421, -1	-1, 23
...	...

Notice how the "-1" seems to only ever appear on the right in ex2 and only ever on the left in "spiky"

[jercgreen6]

@sgbaird
The results are the same:

We remove any of the face labels that feature a -1 before the smoothing occurs:

f = find( any( fl==-1 | fl==0, 2 ) );
fl( f, : ) = [];
tri( f, : ) = [];

[sgbaird]

I was really hoping that would fix it. Taking a second look, I guess that also makes sense that those rows are removed, so the order probably wouldn't matter.

My next best guess is to manipulate SharedTriList.txt. There may be something strange going on with the order. Notice the difference between the first few lines (same pattern across the board as far as I can tell):

ex2	spiky
0, 1, 2	0, 2, 1
1, 2, 3	1, 2, 3
0, 2, 4	0, 4, 2
...	...

This could be messing with some notion of which direction the plane normals are pointing within the code.

You can flip columns 2 and 3 via tri(:, [2, 3]) = tri(:, [3, 2])

which would cause the updated script to be:

folder = fullfile('..', '..', 'examples', 'ex2');

xdat = load( fullfile(folder, 'SharedVertexList.txt') ).';

tri = 1 + load( fullfile(folder, 'SharedTriList.txt') );
tri(:, [2, 3]) = tri(:, [3, 2]);

fl = load( fullfile(folder, 'FaceLabels.txt') );
fl = fliplr(fl); %probably unnecessary, but might as well leave it to make the versions consistent

f = find( any( fl==-1 | fl==0, 2 ) );
fl( f, : ) = [];
tri( f, : ) = [];

ntype = load( fullfile(folder, 'NodeType.txt') );

whos

xsmooth = HierarchicalSmoothMatlab( tri, xdat, fl, ntype );

t = num2cell(xsmooth, 2);
trisurf(tri, t{:}, 'EdgeColor', 'none')
axis equal

Still same issue?

[jercgreen6]

@sgbaird
I get the same plot as before. Additionally, I tried sorting the tri data tri = sort(tri,2);

I also tried multiplying the first and third columns of xdat to match the signs of the example file. This rotates the cube as expected while still maintaining spikes.

xdat(1,:) = -xdat(1,:);
xdat(3,:) = -xdat(3,:);

[sgbaird]

@jercgreen6

Darn..

Here are some other observations I made:

NodeType.txt

It's maybe worth noting that there are 346730 nodes in ex2 and 13886 nodes in ex1.
The unique NodeTypes are the same between the two:

2
3
4
12
13
14

The histograms of the various NodeTypes are also similar:

But there's not as many of NodeType "14" (Quadruple Point on the outer surface) in ex2 as there are in "spiky", and more of "12".

SharedVertexList.txt

In "spiky" the values are all positive. In ex2, some values are negative. I doubt this would make a difference. Just has to do with the center of the microstructure. You tried changing this and had the same issue though.

Other Comments

Try with a similar number of grains/resolution as in ex2. Do you get the same behavior?

Here's a look at the input microstructures:

"spiky" obviously has fewer grains but the resolution is probably about the same and I don't notice anything wildly different qualitatively.

@siddharth-maddali, looks like the code works fine for the examples, but not for the data @jercgreen6 is working with.

[sgbaird]

@siddharth-maddali What's still particularly disconcerting is that I verified the pure MATLAB code seems to work just fine on @jercgreen6's data:

His project involves very large microstructures which would probably make using the pure MATLAB code computationally intractable.

[sgbaird]

Ran it on ex2 as well, and same, works just fine using pure MATLAB code:

[sgbaird]

In other words, since formatting of the input files doesn't seem to be an issue (based on the troubleshooting here), this suggests a bug with the C++ code implementation.

Just for an easy download for myself: spiky.zip

[siddharth-maddali]

@jercgreen6 @sgbaird Can you try out the following and see what happens:
Create a new data set by isolating the set of grain boundaries corresponding to a single grain (that shows up as spiky), and try smoothing that one grain. Do the spikes still show up? You might have to relabel a few things, can be easily done in Matlab.

[sgbaird]

@siddharth-maddali, good suggestion. Thank you for getting back on this. I got my supercomputer account renewed and was able to run the MEX function. It still appears spiky, whether smoothing the full microstucture or a single grain.

Select Spiky Grain

Code

I selected grain #10 and renumbered the tri3 matrix:

i = 10;
ids3 = any(fl2==i, 2);
tri3 = tri2(ids3, :);
fl3 = fl2(ids3,:);
tri3u = unique(tri3);
xdat3 = xdat2(:, tri3u);
ntype3 = ntype2(tri3u,:);

nmax = max(tri2,[],'all');
npts = length(tri3u);
for i = 1:npts
	tri3(tri3 == tri3u(i)) = i+nmax;
end
tri3 = tri3 - nmax;

Data

Here's the data for grain 10: subgrain.mat direct download, with the 3 dropped from the variable names (i.e. xdat, tri, etc.). Triangle 155 tri(155,:) is an example of a jutting triangle, which consists of points 65, 100, and 75 from xdat.

Full Microstructure Smoothing

(left) Original, (right) Pure MATLAB:

(left) Original (right) MEX:

Single Grain

The following is after condensing to a single grain, and then smoothing just the one grain:

Pure MATLAB:

xsmooth3 = HierarchicalSmooth( xdat3, tri3, fl3, ntype3 );

MEX:

xsmooth3 = HierarchicalSmoothMatlab(tri3, xdat3, fl3, ntype3 );

[siddharth-maddali]

@sgbaird The first thing I noticed with the single-grain data set that you shared, is that running the C++ code gives me a warning that about 1/3 of the nodes weren't smoothed. This is concerning, because it means the algorithm couldn't resolve much of the nodes in that dataset, into a well-behaved hierarchy. By this I mean interior (albeit pixelated) surface nodes, bound by triple line nodes on the outer rim, in turn bound by quad points. This is usually a result of faulty labeling of node types (faulty with respect to the algorithm!).

Here is an image of the migration of 'smoothed' nodes from their original position, as you can see very few of them have moved.
The black markers are interior GB points, the red are triple lines and the cyan are quad points addording to your single-grain data set. Magenta arrows show the migration of the individual points due to the smoothing. As you can see, the majority of points don't move, because there is no clear hierarchical relation between the nodes. The C++ algorithm simply ignores the nodes that it cannot put into the hierarchy. I'm guessing you're seeing spikes because it found some "orphaned" interior GB nodes (black markers), in that they had no constraining triple line nodes, and were free to float around. The fact that Matlab seems to be working is confusing to me, I need to look into that one. I suspect that back in the day I implemented a much less stricter rejection policy for the nodes. In the mean time, if you can manually relabel the misbehaving node types to be topologically well-behaved, it will work.

Is it possible to work with a microstructure with finer discretization? This would resolve a lot of these topological problems you're seeing. I realize this is hard, seeing as real-world data is what it is...!

[jercgreen6]

@sgbaird @siddharth-maddali

Thank you for your help and suggestions!

It would be difficult to make the mesh more fine, but I might be able to increase the resolution slightly. The application I am working on requires a collection of twenty 100,000 grain microstructures. Currently, running a microstructure through optimized code can take up to a week on a supercomputer. Therefore, due to time constraints I don't have a lot of wiggle room on resolution.

I did find that there are instances were Dream3D makes microstructures which contain interior grains. This forms an interface between only two grains so I could see how that would potentially be a problem with the labeling system. I was able to remove the interior grains from the data set by adding the following to the code, but I didn't see any change in the smoothing results (the microstructure was still spiky).

Interior Grain (Before Being Removed)

folder = fullfile('..', '..', 'examples', 'ex2');

xdat = load( fullfile(folder, 'SharedVertexList.txt') ).';

tri = 1 + load( fullfile(folder, 'SharedTriList.txt') );
tri(:, [2, 3]) = tri(:, [3, 2]);

fl = load( fullfile(folder, 'FaceLabels.txt') );
fl = fliplr(fl); %probably unnecessary, but might as well leave it to make the versions consistent

f = find( any( fl==-1 | fl==0, 2 ) );
fl( f, : ) = [];
tri( f, : ) = [];
ntype = load( fullfile(folder, 'NodeType.txt') );


%% Remove Interior Grains     
nUniqueFaceLabels = unique( fl, 'rows' );
        
neighbors = NaN(max(fl,[],'all'), mode(nUniqueFaceLabels,'all'));
for j = 1:size(nUniqueFaceLabels,1)

      % Row of "neighbors" is the GrainID. Columns keep track of which grains the grain borders. Padded with NaNs. 

      neighbors(nUniqueFaceLabels(j,1),find(isnan(neighbors(nUniqueFaceLabels(j,1),:)),1)) = nUniqueFaceLabels(j,2);
      neighbors(nUniqueFaceLabels(j,2),find(isnan(neighbors(nUniqueFaceLabels(j,2),:)),1)) = nUniqueFaceLabels(j,1);

end
        
for i = 1:size(neighbors,1)
      % It is an interior grain if it only borders 1 other grain
      if sum(~isnan(neighbors(i,:)),2) == 1
                
             % Delete Shared faces and labels                
              f = any(fl == i,2);
              fl(f,:) = [];
              tri(f,:) = [];
         
       end
end

whos

xsmooth = HierarchicalSmoothMatlab( tri, xdat, fl, ntype );

t = num2cell(xsmooth, 2);
trisurf(tri, t{:}, 'EdgeColor', 'none')
axis equal

I will take a look at editing node types, using other code that I have.

[siddharth-maddali]

Thanks @sgbaird @jercgreen6 for pointing out the potential bug, I will definitely look into it!

[sgbaird]

@siddharth-maddali thanks for the reply! I know this is a rather old project and things can get very busy. Appreciate you taking the time.

@jercgreen6 what about for the test microstructure you have? It would be good to verify that increasing the resolution is one way to resolve the behavior just to make sure nothing else strange is going on (though probably not the solution for your use case). Is that feasible with the pipeline you have?

[siddharth-maddali]

@sgbaird No problem, happy to help!

[siddharth-maddali]

You could revert to Laplacian smoothing for the misbehaving nodes, it doesn't discriminate between node types. But of course, you'll lose the benefit of faithful smoothing...