[kandelabr]

Yes, I check my email like once per day so guess it was a coincidence


About symmetry, classy can't know about that so you have to tell it. Would be quite a thing if it could guess that for an arbitrary sketch!

[Mars409]

Now I would like to try 3D modelling a shrouded axial fan, beginning with change of AoA up to 20 degrees from hub to tip, with a little clearance to the shroud.

Do you have an existing example to show how it can be done?

Am thinking of modifying the exisitng 2D airfoil tutorial, using 0 and 20 degrees AoA for the first and second sketches, respectively, and somehow (need to figure out) create a "cap" to model the tip-to-shroud air gap.

One idea for the "cap" is to
(1) copy the 2nd sketch post optimization to a 3rd sketch, chop up the region within the 20degree tilted airfoil curve by propagating from the pressure- and suction-sides -- have to ensure these sides touch the same number of blocks and the corresponding blocks from opposite sides have identical number of grading cuts -- and from the leading edge to the trailing edge;
(2) copy the 3rd sketch to a 4th sketch and translate it radially outwards a little to touch the shroud wall; and
(3) loft from the 3rd to the 4th sketch, employing multigrading in the radial (i.e. 3rd) axis to apply fine grading adjacent to the 3rd sketch.

One serious drawback of this cap is the fine gradings near the 3rd sketch will be touching coarse gradings from the air region before this cap.

A better way would be to create a rounded dome to cap the tip end and somehow make chops that wrap from the leading edge to the trailing edge. Haven't thought of how this can be done. Maybe you already have done this?

[kandelabr]

I've had this in mind for quite a while but never had the time.


I don't imagine exactly how this "cap" should look like, do you have a sketch or a reference mesh?


If I was to do this, I'd first upgrade the Sketch so that it can "wrap" onto a cylinder at that radius, then just loft between sketches. Use the same clamps as in the tutorial (syntax is a bit different, check the latest examples) but optimize in 3D as this is safer.


What about AMI/NCC or other kinds of mesh separation/joining? It could simplify this task greatly (again, I can't tell because I don't know what exactly you look for...)

[Mars409]

After my last post, I thought about the coarse-fine axis-2 grading split across the 2nd sketch and realized it could be avoided by applying multiple-grading to the loft btwn the 1st and 2nd sketches to get progressively finer and finer cell thickness towards the 2nd sketch.

It's all still in the mind.

Next I try to see in my mind how to make a smooth "dome" to cap off the tip and how to make the boundary layer gradings in the previous lofts that have edges drawn on the airfoil of the 2nd sketch wrap over this "dome".

I will use cyclic boundary (instead of NCC) to model just one blade out of N blades of the axial fan. So this places an optimization constraint that points 1 and 8 must share the same x coordinate value and so do points 2 and 7 (How do you apply a constraint like this in Classy_block?).

But I don't see how this (cyclic boundary) makes the job of creating a script for the cap region any easier.

I don't have a visualization in the mind of how your hint ("upgrade the Sketch so that it can "wrap" onto a cylinder at that radius") works. Is it the 2nd sketch that is to be upgraded to wrap onto the cylinder as the end-cap? Or is it the 3rd sketch? Perhaps you can elaborate further.

2nd thoughts about wrapping tip end with boundary layer
===========================================
A day after, on further thoughts, it turns out not necessary to wrap the tip end of the fan blade in boundary layer, because the blade thickness is small and so any viscous forces on the fan and their torque is negligible. So it would be alright to let the boundary layers in the 2nd sketch to loft all the way to the shroud. IOW, the 3rd (at tip end) and the 4th (at shroud wall) sketches will both have the boundary layers like in the 1st & 2nd sketches.

The 3rd & 4th sketches, however, will have blocks filling out the inside of the airfoil curve. These fill-in blocks require no new points -- their top and bottom sides are on the airfoil curve's suction and pressure sides, respectively.

The block at the trailing edge of the airfoil in the 2nd sketch thus in the 3rd sketch continues left into the airfoil's inside across these fill-in blocks. At the left end of the airfoil, i.e. the leading edge, however, to continue, a new block left of the airfoil leading edge needs to be inserted between the existing blocks 0 & 6.

Let's number this newly inserted block 10. It stretches from the airfoil leading edge to the left side of the 2D domain. Two new points, 18 & 19, are inserted between points 12 & 13 and between points 0 & 8, respectively. The newly inserted leading edge block 10 thus consists of points {0, 12, 18, 19}.

With newly inserted point 18 (on the airfoil leading edge), we can form a fill-in block that consists of points {18, 13, 11, 12}. (Basically, there needs to be an integer multiple of pair of points on the airfoil curve in order to completely fill out the inside of the airfoil curve with blocks.)

Then, the next two fill-in blocks consists of points {13, 14, 10, 11} and points {10, 14, 15, 9}.

Let's number these 3 fill-in blocks 11, 12, 13. Together, the newly inserted leading edge block 10, these 3 fill-in blocks 11, 12, 13, and the existing trailing edge blocks 3, 8 form a row of blocks and so for axis 2 chop specifying the new block 10 suffices to determine the left-to-right (with respect to the sketch) chopping for these blocks.

And top-to-bottom (wrt to the sketch) chopping in axis 1 for row of blocks {0, 11, 6}, row {1, 12, 5}, and row {7, 2, 13, 4, 9}.

Now, in view of the chopping axes, the quads for these new blocks are ordered as below:
10 : {18, 12, 0, 19 },
11 : {11, 12, 18, 13 },
12 : {10, 11, 13, 14 }, and
13 : {9, 10, 14, 15 },
while old block 6 becomes:
6 : {13, 18, 19, 8}

Now, obviously we cannot have both Sketches 2 & 3 sit at the same z coordinate. So Sketch 2 is no longer used for lofting but Sketch 1 directly lofts to Sketch 3 with the proviso that there be no loft for faces within the airfoil curve.

But there's also the newly inserted leading edge block 10 that needs to back into Sketch 1, which means Sketch 1 too needs to have the new block 10. Let's call this revised Sketch 1A.

And so Sketch 1A lofts to Sketch 3, except the blocks within the airfoil are excluded from the lofting (Does the existing lofting supports this? Answer: Yes. Each loft is constructed for a quad of points, so just iterate thru the quads we want and skip those within the airfoil curve. Cf function get_loft() in airfoil.py.).

And Sketch 3 is copied as Sketch 4 and displaced further out a short distance, say 1mm, in the z-direction.

And Sketch 3 lofts to Sketch 4.

In all of Sketch 1A, 3 and 4, each sketch needs to constrain the outermost points on the suction side to overlap exactly the corresponding points on the pressure side under the requisite rotation about the fan rotation axis in order to work with cyclic boundary. Is there such an optimization constraint existing in Classy_blocks?

[kandelabr]

Yes, there's RotationLink (https://github.com/damogranlabs/clas...e/links.py#L44) just for what you need (or, I need for pump impeller generation).


Your plan could work and it would be simple to try it out with the latest version. Once you have a sketch, you can simply do

shape = cb.ExtrudedShape(sketch, amount) 
# or
cb.LoftedShape(sketch_1, sketch_2).

As long as you're just developing blocking, you'll need coarse cell counts and a simple chopping method that won't get in the way. You can use

grader = cb.FixedCountGrader(mesh, cell_count)
grader.grade()
# or
grader = cb.SimpleGrader(mesh, cell_size)
grader.grade()


# just before
mesh.write(...)

Some voodoo will be required to properly block the innermost dome (if I understand you correctly). You can deal with that later once you have your blade blocking figured out.


Let me know how it goes.

[Mars409]

I am looking for an optimization constraint that forces identical x coordinates btwn points 2 & 7, b/c the 1-2-3 boundary transform onto the 6-7-8 boundary under cyclic boundary, the rotation being around the fan axis which runs left-right some radius behind the page. So, point vertex 2 can slide left-right along the 1-2-3 boundary, and vertex 7 will have to follow. In other words, vertex 7 does not have any degree of freedom.

RotationLink, which keeps the follower at the same radius given a rotation origin and a normal to the plane of rotation, on the other hand, rotates the follower by the same angle that the leader has been rotated. It's different.

SymmetryLink, which makes the follower a symmetry image of the leader across a given plane, would work if my sketch is symmetrical about the x axis, which isn't true in most cases.

So I am going to look at the optimization codes to see how the LinkBase class is used. My guess is that each instantiation of the LinkBase class for a follower point removes the follower point from a list of points to optimization position of, and the transform() method is called to compute the follower position at each change of the leader position. If so, then it suffices to create a link class derived from LinkBase class for which the transform() method merely changes the follower's x coordinate to the leader's.

After that, the plan seems complete, and time seems ripe to go for a try.

There's just one more problem. Strictly speaking, each sketch lies on the surface of an imaginary cylinder whole centerline is the centerline of the axial fan. And lofting should be radiating from Sketch 1A radially outwards to Sketch 3. So, in fact, Sketch 3 is proportionately larger than Sketch 1A by the ratio of the radial distances from the fan axis. In short, the cylindrical coordinate system is more suitable for the axes than the cartesian. Are there classes of Classy_blocks that support sketch, loft, and mesh assembly in cylindrical coordinate system?

Along this line, to support cyclic boundary, I would also like to have all points on the lower boundary of each sketch to translate vertically from the upper boundary of the same sketch by the same rotation angle for all lower sketch boundary points of all sketches. Any existing method to support this?

[kandelabr]

No, there is no such Link implemented but it will be easy to make one, as you have found out.


You "last" problem is the issue I talked about before, that is, "wrapping" a sketch on a cylinder plane. I haven't thought about this yet. I guess the simplest way would be to create a polar sketch that would transform to cartesian after everything has been done. That being said, it will require some handling for transformations and whatnot... Your scaling a.k.a. "radius increase" would then be solved automatically. I'll think about it.


I had problems with cyclic BCs and blockMesh in the past and never really solved them because I had to create entire geometry. Just now I'm working on a simplified case with cyclics and I'll have to sort it out. Will update the package once I have this figured out and let you know.

[Mars409]

I think I will just project Sketch 1A, 3, and 4 onto the cylinders, all centered around the fan rotation axis, at the respect radii and loft between them, instead of lofting directly between those flat sketches.

And to enforce cyclic boundary between the upper and lower boundaries of each projected sketch, points on the upper boundaries of the sketches post projection will be transformed to the respective lower boundaries by rotation.

But since there's no rotation transformation in blockMesh post projection, the rotational transformation to properly position the lower boundary points should be done first, followed by projection to the respective cylinders.

This, however, runs a risk of loss of accuracy b/c the resulting vertices may not map exactly onto the upper boundary vertices post projection, especially in case the blade tip is circumferentially swept forward significantly, e.g. by 80 degrees around the fan axis.

But on 2nd thought this may not be a problem if Sketch 3 is rotated by 80 degrees around the fan axis before projection onto the cylinder.

But no, this will still not work where the sketch is non-rectangular, e.g. a rhombus shape with slanted upper & lower boundaries. For example, in the rhombus case, and left boundary being lower and the right boundary being higher in order to follow the pitch (angle the chord makes with the fan axis) of the airfoil, the sketch's upper left and lower left corners will be at different distances from the cylinder and, as a result, their radial projections onto the cylinder may not rotationally transform into each other.

And the correct positions of the lower boundary points need to be available to the optimizer class anyway.

===== Conclusion ===

So, the safer way to guarantee cyclic boundary btwn the upper and lower boundaries of the sketch post projection is still to do the rotational transformation from the former to the latter post projection.

That leaves us to do the projection of points in Sketch 1A, 3 and 4 onto the cylinders at the respective radii within Python, followed by the rotational transformation of the upper boundary points around the fan axis to give the correct positions for the corresponding lower boundary points.

A separate axial rotation is then applied to the post projection sketch to the final position to effect the circumferential sweep.

What's left now is to figure out how to incorporate these changes into Classy_blocks. Am beginning to think about this.

[kandelabr]

Projections should be exact if you define a searchableCylinder. They definitely will not be if you project to an STL file because each blade will get a different set of triangles. That's not a problem though, because you can define extra planes through specified points and then project to both surfaces, plane and cylinder; this will ensure the curves will be the same.


Or to avoid projections at all and create a spline edge with pre-calculated points - intersection between a cylinder and blade. I'm working on something similar right now...


Based on my experience it's better to keep away from projections if possible.

[Mars409]

Here is my Sketch3 at the tip of the blade. All block and point numberings follow your Tutorial#6 airfoil sketch except new blocks 10-12 are added (3 instead of 1) on the left and new blocks 13-16 are created within the airfoil curve. The left boundary of the sketch is a straight line (in the outlet plane of rotation) instead of semicircle.

Without considering projection onto cylinders nor change of pitch (the angle between the chord of the airfoil and the left-right axis of the sketch) for the time being, I will simply extrude into the page a distance to the hub for all blocks BUT those within the airfoil curve, and extrude ALL blocks out of the page about 1/4*CHORD. Then I will optimize, backport, discard the out-of-page extruded lofts, and redo out-of-page lofting, this time to precisely loft a small distance, around 1mm, with no optimization to follow.

The first out-of-page lofts were there to nicely optimize the gradings in the faces 13-16. Their thickness 1/4*CHORD is chosen to avoid aspect ratio issue arising from very thin lofts during optimization. On the optimization is done, these sacrificial out-of-page lofts are discarded in favor of the real final 1mm-thick beyond-tip lofts.

What objects are backported that I can use as argument to build lofts from? My first guess would be the vertices (i.e. points with indices) but need to inspect the Python code to be sure.

Once I get pass this, then I will play with in-code projection to imaginary cylinders at the 3 different radii (hub, tip, and shroud), followed by rotating points on upper sketch boundaries around the fan axis (extending left-right a radius behind the sketch) to replace corresponding points on the respective lower sketch boundaries to honor cyclic boundary.

There's also the linear scaling of the sketch top-bottom dimension in proportion to the radii. I am not sure about this one. It seems to have a disadvantage of disturbing the geometry of the airfoil curve. Maybe I will just define the y coordinates of the points on the upper sketch boundary in direct proportion to the radius (and the lower sketch boundary points follows under the abovementioned rotation.

[kandelabr]

When you add Shapes/Operations to the mesh, nothing is really done. Before running MeshOptimizer, the mesh needs to be assembled, meaning operations' points are transformed into vertices and they are given indexes. Optimization moves those vertices (changing Blocks) but original operations remain the same. Mesh.backport() moves original operations so that they correspond to blocks after optimization.


If you use ShapeOptimizer or SketchOptimizer, original operations are moved even before they are added to the mesh.

[Mars409]

I see bottom_face and top_face of Face class are attributes of the operation class, and points is an attribute of the Face class initialized from a Python list of points (PointListType) using Numpy asarray() method so any changes made by the MeshOptimizer to the vertices indeed is reflected in the original points, e.g. the airfoil.py points ndarray created in line 73. So I can just refer to points by index in this points list instead of by way of block and top/bottom face and then calling points() method.

Edit
===
Oh wait, I was wrong. Face class' points attribute is its own attribute, not reference:

self.points = [Point(p) for p in points]

and so, backport() will only update the points in the faces of the original operations, but not the original points list such as the ndarray created in line 73 of airfoil.py. As such, to access the optimized points, after mesh.clear(), the only way is thru the original operations.

SEO suggestion
=============
Classy_blocks didn't show up when I initially searched for "Openfoam structurered meshing", but only much later when I looked for "airfoil mesh."

You might want to add "structured mesh" to the list of search terms on your Github pages and that of damogranlabs.com.

[kandelabr]

You are right, points are an attribute of Faces and Operation's points are just a read-only shortcut to faces' data. You can't modify operation's points directly. Furthermore, if you would want to change points in a Shape, Stack or Object, you'd have to manually find all other coincident faces and change their points too or they will remain unchanged. That would be a mess, of course, so you can either use an optimizer of some kind or make a Grid from your entities:

from classy_blocks.optimize.grid import QuadGrid, HexGrid


hex_grid = HexGrid.from_elements(...)
# hex_grid.positions is a property and you can't change that directly


# change a point like that:
hex_grid.update(index, position)


# similar for a flat grid:
quad_grid = QuadGrid.from_sketch(...)

Thanks for the suggestion! I'll update keywords and tags.

[Mars409]

Am looking at Face.add_edge() method.

Is it true that add_edge() must be called with the same edge data specification for all faces that share the same edge? Or just doing it for any one of the abutting faces will do?

What if somehow the same edge is given different edge specification in abutting faces from abutting operations?

I guess during mesh.assemble() and list of edges is gathered from all the operations that have been added to the mesh.

So, here, will the doubly specified edge just quietly adopt the last defined edge data for the edge without issuing warning of conflict of specifications?

[kandelabr]

The first edge will be taken into account. When assembling the mesh, the script first searches for an existing edge and adds a new one when nothing was found.


https://github.com/damogranlabs/clas...ge_list.py#L26

[Mars409]

Thanks.

The python script is almost there.

Testing & debugging is next. Could be challenging.

[Mars409]

Have been going thru the mesh and optimizer codes.

Some names confuse me everytime.

https://github.com/damogranlabs/clas...ks/mesh.py#L67 :

new_vertices = self.vertex_list.add(point, list(patches))
vertices.append(new_vertices)

Here, new_vertices really is just a single vertex, not a list of multiple vertices. vertice, in the singular, seems more appropriate.

https://github.com/damogranlabs/clas...s/mesh.py#L191

vertices = [vertex.position for vertex in block.vertices]
op.bottom_face.update(vertices[:4])
op.top_face.update(vertices[4:])

Here, vertices on the LHS of line 191 seem more appropriately names points, b/c vertex is position plus index, whereas point is just position.

====

Changing point positions by way of user's operations
========================================
Looking at the grid codes, it appears not quite handy for the user to identify which vertex is which target point in the target operation to be edited. In my present fan blade model, I need to edit the points on the last sketch to place them on the shroud surface after optimizing and backport. The grid doesn't provide the hooks to find the vertices that correspond to those points on the shroud surface, so despite it providing update() as the hook I don't have the vertex index to give it.

Fortunately, after backport, the vertices positions are updated in the points in user's operations. They can be read by way of operation.bottom/top_face.point_array(), and can be edited by way of operation.bottom/top_face.update(points).

The drawback is the same point must be consistently updated with the same position value in every one of the operations that it is part of. Here is where the grid.update() has an advantage, b/c it needs to update the id'ed vertex once.

So, seems this is how I can proceed.

====

Changing vertex position by way of grid
==============================
Well, that may not be the last word, b/c the vertex ids can be read from the mesh by way of mesh.blocks[i].vertices[j]:

vertex_id = mesh.blocks[i].vertices[j]

and here i is known b/c it is the sequence number of the corresponding operation in the sequence of calling mesh.add(operation) and 0 <= j < 4 is known for each of the corners of the operation of each of the bottom/top faces.

To use this route (grid.update()), this method will be called after calling MeshOptimizer(mesh) but before calling mesh.backport().

So, this approach looks more attractive after all.

And, even in case the vertex editing needs to happen after the mesh backport has happened, we still can rebuild the grid from the mesh in order to access grid.update():

grid = HexGrid.from_mesh(mesh)

On further investigation, turns out just accessing the existing grid attribute in the optimizer will do. Just update the target vertices, then call backport() :

optimizer.grid.update(vertex_id, new_position)
optimizer.grid.backport()

Even in case no optimizer has been created to actually perform optimization, one can be created just for this purpose and no more :

optimizer = cb.MeshOptimizer(mesh)

[kandelabr]

Yes, Grid was not meant to be used directly by the user so there are no 'convience' methods. But if it turns out it is more useful than its current scope, those methods should really be added.

[Mars409]

Still working on the fan blade Python codes. One question arises: If grading is not completely specified due to omission in the operation.chop() calls, which class issues a warning?

[kandelabr]

Precisely here: https://github.com/damogranlabs/clas...ck_list.py#L67