AssertionError

[WoolenWolfff]

Hello!
If I run grain_neighborhoods.py once from the IDE by simply pressing "run" - everything works fine, a new image is generated each time. But if I want, say, to run it several times in a row, like this:

def main(num):
# here is the code from rain_neighborhoods.py

for i in range(5):
    main(i)

As a result, I do not get images with the names "trimesh0.png", "trimesh1.png", ..., but get an AssertionError. The same thing happens if you try to run a cell several times in jupyter:

---

AssertionError Traceback (most recent call last)
Cell In[2], line 79
76 repl_seeds.append(ring_seed)
78 #Create polygon and triangle meshes
---> 79 pmesh = msp.meshing.PolyMesh.from_seeds(repl_seeds, domain)
80 phases = [{'material_type': 'solid'} for i in range(4)]
81 phases[0]['material_type'] = 'matrix'

File
\microstructpy\src\microstructpy\meshing\polymesh.py:639, in PolyMesh.from_seeds(cls, seedlist, domain, edge_opt, n_iter, verbose)
636 removing_pts = False
638 missing_seeds = set(range(len(seedlist))) - set(bkdwn2seed)
--> 639 assert not missing_seeds, str(missing_seeds)
640 # compute voronoi diagram
641 voro_fun = {2: pyvoro.compute_2d_voronoi,
642 3: pyvoro.compute_voronoi}[n_dim]

AssertionError: {2, 11, 1038, ....}

---

If I understand correctly, the problem is that "missing_seeds" is not empty. I just don't understand how to fix it.
In the future, I will need to generate many different structures with different parameters. I would like to automate this. But until this problem is solved, I'm afraid it will not be possible.

Please help me solve this problem!

[kip-hart]

Hi @WoolenWolfff , this error occurs because sometimes the Laguerre-Voronoi tessellation does not contain every seed particle. This sort of thing creates issues when trying to trace from seed to polygon to mesh element. That's why I threw the AssertionError in there, to avoid that issue.

In terms of getting your script to run, I would try to embed the for loop within main, rather than call main 5 times. I don't know what gets saved locally vs restarted new on each call.

[WoolenWolfff]

You write "I would try to embed the for loop within main, rather than call main 5 times". The thing is, I tried different ways, the point is that I can't get a normal result without an error if I run the code several times. Let's say even if I run it from a jupyter notebook cell, the code only runs once, then I get an AssertionError.

Maybe there is some way to completely clear the seeds before a new generation?

[kip-hart]

Could you post the file you're running so I can attempt to reproduce?

[WoolenWolfff]

I just copy the code from the grain_neighborhoods.py example into a cell and try to run it more than once. The first time it works, but then I get an error. I can't attach the file, but it literally has this code in one cell:

from __future__ import division

import os

import numpy as np
import scipy.integrate
import scipy.stats
from matplotlib import pyplot as plt

import microstructpy as msp

# Define the domain
domain = msp.geometry.Square(corner=(0, 0), side_length=10)

# Define the material phases
a_dist = scipy.stats.lognorm(s=1, scale=0.1)
matrix_phase = {'fraction': 1,
                'material_type': 'matrix',
                'shape': 'circle',
                'area': a_dist}

neighborhood_phase = {'fraction': 1,
                      'material_type': 'solid',
                      'shape': 'ellipse',
                      'a': 1.5,
                      'b': 0.6,
                      'angle_deg': scipy.stats.uniform(0, 360)}

phases = [matrix_phase, neighborhood_phase]

# Create the seed list
seeds = msp.seeding.SeedList.from_info(phases, domain.area)
seeds.position(domain)

# Replace the neighborhood phase with materials
a = neighborhood_phase['a']
b = neighborhood_phase['b']
r = b / 3
n = 16

t_perim = np.linspace(0, 2 * np.pi, 201)
x_perim = (a - r) * np.cos(t_perim)
y_perim = (b - r) * np.sin(t_perim)
dx = np.insert(np.diff(x_perim), 0, 0)
dy = np.insert(np.diff(y_perim), 0, 0)
ds = np.sqrt(dx * dx + dy * dy)
arc_len = scipy.integrate.cumtrapz(ds, x=t_perim, initial=0)
eq_spaced = arc_len[-1] * np.arange(n) / n
x_pts = np.interp(eq_spaced, arc_len, x_perim)
y_pts = np.interp(eq_spaced, arc_len, y_perim)

repl_seeds = msp.seeding.SeedList()
geom = {'a': a - 2 * r, 'b': b - 2 * r}
for sn, seed in enumerate(seeds):
    if seed.phase == 0:
        repl_seeds.append(seed)
    else:
        center = seed.position
        theta = seed.geometry.angle_rad

        geom['angle_rad'] = theta
        geom['center'] = center
        core_seed = msp.seeding.Seed.factory('ellipse', phase=3,
                                             position=seed.position, **geom)
        repl_seeds.append(core_seed)

        x_ring = center[0] + x_pts * np.cos(theta) - y_pts * np.sin(theta)
        y_ring = center[1] + x_pts * np.sin(theta) + y_pts * np.cos(theta)
        for i in range(n):
            phase = 1 + (i % 2)
            center = (x_ring[i], y_ring[i])
            ring_geom = {'center': center, 'r': r}
            ring_seed = msp.seeding.Seed.factory('circle', position=center,
                                                 phase=phase, **ring_geom)
            if domain.within(center):
                repl_seeds.append(ring_seed)

# Create polygon and triangle meshes
pmesh = msp.meshing.PolyMesh.from_seeds(repl_seeds, domain)
phases = [{'material_type': 'solid'} for i in range(4)]
phases[0]['material_type'] = 'matrix'
tmesh = msp.meshing.TriMesh.from_polymesh(pmesh, phases, min_angle=20,
                                          max_volume=0.1)

# Plot triangle mesh
colors = ['C' + str(repl_seeds[att].phase) for att in tmesh.element_attributes]
tmesh.plot(facecolors=colors, edgecolors='k', linewidth=0.2)

plt.axis('square')
plt.xlim(domain.limits[0])
plt.ylim(domain.limits[1])