CHLNDDEV / OceanMesh2D

A two-dimensional triangular mesh generator with pre- and post-processing utilities written in pure MATLAB (no toolboxes required) designed specifically to build models that solve shallow-water equations or wave equations in a coastal environment (ADCIRC, FVCOM, WaveWatch3, SWAN, SCHISM, Telemac, etc.).
https://github.com/sponsors/krober10nd
GNU General Public License v3.0
179 stars 65 forks source link

fix to carry over weirs #206

Closed krober10nd closed 3 years ago

krober10nd commented 3 years ago
WPringle commented 3 years ago

Are the original arrays guaranteed to be row vectors? Why we have a discrepancy here?

krober10nd commented 3 years ago

with all the manipulation in the plus perhaps they become rows somewhere there.

WPringle commented 3 years ago

Not sure about that. What is the canonical format that we are using? (row or column vector)

krober10nd commented 3 years ago

everything should be column vectors.

WPringle commented 3 years ago

Yeah so for some reason you are saying that originally "obj1" has them as a row vector. Was it like this before the call to function? Or did it get changed during plus function?

Need to see where this happened and prevent it there instead.

krober10nd commented 3 years ago

Not sure. I don't see why it's a big deal. Just make sure they're columns by doing (:) when carrying over.

WPringle commented 3 years ago

Yeah but the transpose operator may make the shape invalid for another test case. Perhaps it was just wrong in your original mesh

krober10nd commented 3 years ago

No, it wasn't wrong it my original mesh. It should be column vectors so this PR forces that to be true. What other operations are your referring to that would make it so-called "invalid"?

WPringle commented 3 years ago

You are using transpose operator to force a row vector to a column vector. But it may already be a column vector. So the PR doesn't necessarily force it to be a column vector (it forces a transpose).

krober10nd commented 3 years ago

That's true. So it needs to use a (:).

krober10nd commented 3 years ago

Try to plus m4 into m2

m3 = plus(m4,m2,'arb',{'ds',0});

both m4 and m2 have weirs. Specifically m2 the parent has 6 weirs, and m4, the child, has one.

These meshes were created using all the available tools here in om.

edit: I'll send you the file via email. it's too big to upload

krober10nd commented 3 years ago

For example in the case I sent to you via email, the logical array jj is 385 entries long

   obj.bd.ibtype = [obj.bd.ibtype ; obj1.bd.ibtype(jj)];

and obj.bd.ibtype is a scalar. When indexing into obj.bd.ibtype(jj) it produces a row array. This obviously does not concatenate.

Perhaps this case was not previously encountered in our works.

WPringle commented 3 years ago

Keith, so it seems that make_bc is making the ibtype and nvell arrays on the fly so they are populated as row arrays, i.e. dimensions of [1 x NBOU]. This does happen to be consistent with nbvv having dimensions of [max(nvell) x NBOU].

So the solution is to change the semi-colon when concatenating to just a blank space or a comma rather than transposing.

It so happens that the shape of ibtype and nvell doesn't matter for writing out the mesh. But when reading in a mesh from fort.14 it was assumed to be dimensions [NBOU x 1]. So I made an edit for that too.

krober10nd commented 3 years ago

Yes I had read in the base mesh from a fort14 so perhaps that why this happened. Did the case I sent now work?

krober10nd commented 3 years ago

Okay, I ran the case, everything works. I agree the fix makes sense. Lets merge this.