Extract numerical values from boundaries

[lrvie]

Simply put, how to extract the coordinates and an expression evaluated at a given boundary directly from the model attribute?

More specifically, consider a given boundary of a 2D-axisymmetric component. I want to extract its coordinates and also an expression, say (r,z) and comp1.T.

I haven't found a direct MPh function to do it, but it seems to be possible via the .java attribute. Unfortunately, I couldn't find the answer either.

From comsol's programming reference manual, it seems that the numerical() attribute could do the trick, e.g.: model.result().numerical(<ftag>).getData();, though it's not clear to me how to reference the desired boundary and expression.

Does anyone happen to know a way to approach this?

[john-hen]

If you can't figure this out, consider attaching a Comsol model to your question, for demonstration purposes. Or refer to the demo model from the tutorial. Yes, this is definitely possible if you go through the Comsol API, via the .java attribute. And there might be an even easier way using MPh's Model.evaluate() if you define some kind of "coupling operator" in your model. Hard to answer without a concrete example.

[lrvie]

Thanks for your answer @john-hen. Consider the attached model. What I'm trying to do is pretty much the same as the 1D plot group, i.e. extract the (r,z) coordinates of a deformed boundary (boundary 4) alongside with an expression evaluated along this boundary (here, the colored velocity spf.U).

If it's possible to use model.evaluate it would be great, though I don't know how to reference the desired boundary.

Otherwise, I was experimenting with a java method within Comsol's Application Builder and it's possible to extract velocities at a given point using the at2 operator

NumericalFeature newEval = model.result().numerical().create("newEval", "EvalGlobal");
newEval.set("expr", "at2(0.010,0.014,spf.U)");
newEval.set("data", "dset1");

double[][] data = newEval.getReal();
debugLog(data);

model.result().numerical().remove("newEval");

which gives

{{5.395839106504864E-4, 3.6587242817664306E-4, 5.766007063265596E-5, 2.6563618985988305E-7, 2.706686382308497E-7}}

But again I don't know how to reference the desired boundary instead of a point. Also, translating this to mph's python script will return a <java array 'double[][]'> variable. Maybe one has to parse it somehow?

Hope it provides more details about my question :)
extract_values.zip

[john-hen]

There are probably a number of ways to do this. For example, the Comsol API also lets us export the data associated with a specific plot. But it seems more straightforward to me to define an edge dataset and get to that data directly. There may be a way to use a "linear extrusion" operator, or some such operator, so that we could just call MPh's Model.evaluate(), though I haven't looked into that much. It would not accept an edge dataset, since it uses Comsol's Eval feature under the hood, which is limited to the most common dataset types, like the ones that contain the actual solution coming out of the solver sequence.

Instead, we can use an Interp feature, which accepts arbitrary datasets and "uses interpolation". Not sure exactly "what" it interpolates, or rather what determines the discretization in the end. But it gets you the data along a given boundary, probably determined by the mesh density there.

Consider this script:

import mph
from numpy import array

client = mph.start()
model = client.load('extract_values.mph')

# Create edge dataset and select boundary number 4.
edge = (model/'datasets').create('Edge2D')
edge.select(4)

# Retrieve the data for each field.
interp = (model/'evaluations').create('Interp')
interp.property('data', edge)
interp.property('solnum', 1)      # Only first time step (inner solution).
interp.property('expr', 'r')
interp.property('unit', 'mm')
r = array(interp.java.getReal()).squeeze()
interp.property('expr', 'z')
interp.property('unit', 'mm')
z = array(interp.java.getReal()).squeeze()
interp.property('expr', 'spf.U')
interp.property('unit', 'mm/s')
U = array(interp.java.getReal()).squeeze()

# Print some stats.
print(f'Retrieved {len(r)} values.')
print(f'r = {r.min()} ... {r.max()} mm')
print(f'z = {z.min()} ... {z.max()} mm')
print(f'U = {U.min()} ... {U.max()} mm/s')

Which outputs:

❯ python eval_edge.py
Retrieved 175 values.
r = 10.0 ... 20.0 mm
z = 0.0 ... 15.0 mm
U = 0.0 ... 0.8441769135022498 mm/s

[lrvie]

Sorry for the late reply @john-hen and many thanks for the useful answer.

I've been testing it and it works pretty decently, however I have a follow-up question: If I the Edge 2D dataset already exists in the original model, and there are multiple models, say model1, model2, etc., I assume I don't need to create a new dataset every time. But I haven't been able to use your code in conjunction with existing datasets.

I tried something like:

edge1 = (model1/'datasets/Edge 2D 1')
interp1 = (model1/'evaluations/Interp')
edge.select(4)
interp1.property('data', edge1)
interp1.property('solnum', 2)
interp1.property('expr', 'r')
interp1.property('unit', 'mm')
out = array(interp1.java.getReal()).squeeze()

Which gives the error

---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last) 
Cell 9 line 4
      [2] interp1 = (model1/'evaluations/Interp')
      [3] edge1.select(4)
----> [4] interp1.property('data', edge1)
      [5] interp1.property('solnum', 2)
      [6] interp1.property('expr', 'r')

File ..., in Node.property(self, name, value)
    [394] return get(self.java, name)
    [395] else:
--> [396] self.java.set(name, cast(value))

AttributeError: 'NoneType' object has no attribute 'set'

Can you enlighten me on how to properly refer to existing datasets in this case?

[john-hen]

As I see it, this error occurs because the evaluation node model1/'evaluations/Interp does not exist. If you called interp1.exists() it should return False. Either you never created that node, or it's called something other than Interp. You'll have to inspect your model to find out, either in the Comsol GUI or using mph.tree(). Nodes have the names you give them, or they have the default name Comsol gave them on creation if you did not change it.

[lrvie]

Ok, after inspecting closer using mph.tree(model1) I discovered that the interpolation created had the name Interpolation 1, therefore the command should be interp1 = (model1/'evaluations/Interpolation 1'). Alternately, one can create an interpolation node with a custom name as interp = (model/'evaluations').create('Interp', name='custom name') for easier future reference. Thanks!

[SereinQAQ]

How to export the concentration information on an edge
a = java_model.result().dataset().create("edg1", "Edge2D").selection().set(40)
a = model / "datasets" / "Edge 2D 1"
print(model.evaluate("c", dataset=a))

RuntimeError: Dataset "Edge 2D 1" does not refer to a solution.
@john-hen @lrvie I would appreciate it if you could answer me.

[lrvie]

Without having a look at your model is hard to tell. From the error message it seems that there's no dataset named "Edge 2D 1". I'd also double check if your system requires slash "/" or backslash " \ " to navigate through the model tree. Windows and Linux have different conventions for example. Also try something like
a = model/"datasets/Edge 2D 1" or a = model\"datasets\Edge 2D 1"

[john-hen]

@SereinQAQ
I don't see how your question isn't covered by the answer further up. It explains why Model.evaluate() doesn't work for your use case, but also explains how to get around it.