Hello Bill,
Thank you very much for this! This looks very much like the domain I need to
work with. I will try modifying it to fit the domain and let you know the
results.
Thank you for your time,
Ohi
________________________________
From: overture-bounce@xxxxxxxxxxxxx <overture-bounce@xxxxxxxxxxxxx> on behalf
of Bill Henshaw <wdhenshaw@xxxxxxxxx>
Sent: Monday, June 19, 2017 4:03:30 PM
To: overture@xxxxxxxxxxxxx
Subject: [overture] Re: Matrix has a row with all zero coefficients
Hi Ohi,
Unfortunately your grid is not the greatest. You should basically never
have to use prevent hole cutting
when building a grid. Here is a grid for a back step --- try this and it it
works you can
add more grids to make it look like the domain you want.
Regards
Bill
On Mon, Jun 19, 2017 at 5:21 PM, Ohiremen L Dibua
<odibua@xxxxxxxxxxxx<mailto:odibua@xxxxxxxxxxxx>> wrote:
Hello Bill,
Thank you for your response. I tested the grid on tcm3, and found that the
issue was that I needed to add an extra ghost point. I previously had 1 ghost
point, but I added 2 instead. This allowed my code to run, and the Laplacian to
be solved, but now the solution is one that does not make sense. I have
attached a coarse version of the mesh grid, and the solution to this message.
Although the solution is less accurate, the pattern seen is consistent. There
are 2 regions which have a solution of zero within the domain. These regions
have boundaries that interpolate with the neighbors on every side, but one in
which I impose symmetric BCs (ie du/dn=0). It took me a while, but I believe
the issue is that, as a result of my BC's these 2 regions have an f of zero at
all of the boundaries. This results in a solution of zero within the region.
However, I thought that sharing the boundaries of these regions with the
adjacent grids, and allowing implicit interpolation would prevent this from
occurring. I have spent a lot of time trying different things to fix this, but
have been unable to. I wanted to ask if you had any thoughts on things I should
try in Overture?
Below is the cmd file for generating this grid.
Thank you for your time,
Ohi
create mappings
#Create Domain 1 ###############
rectangle
mappingName
domain1
set corners
-12.5 2.5 0. 1.5
lines
10 10
#100 100
boundary conditions
1 0 1 1
share
0 1 0 0
exit
#Create Domain 4 ###############
rectangle
mappingName
domain4
set corners
2.5 7.5 0 0.5
#-12.5 7.5 0 0.5
lines
10 10
#100 100
boundary conditions
0 0 1 0
share
0 0 0 3
exit
#Create Domain 5 ###############
rectangle
mappingName
domain5
set corners
7.5 12.5 0 1
lines
10 10
#100 100
boundary conditions
0 0 1 1
exit
#Create Domain 3 ###############
spline
enter spline points
5
7. 1.5
7.125 1.169
7.25 1.067
7.375 1.016
7.5 1.
lines
10
#100
mappingName
topSpline
exit
spline
enter spline points
7
7.5 1.
7.5 0.9
7.5 0.85
7.5 0.8
7.5 0.7
7.5 0.6
7.5 0.5
lines
10
#100
mappingName
rightSpline
exit
spline
enter spline points
9
6.5 1.5
6.55 1.188
6.6 1.064
6.7 0.9
6.8 0.789
7. 0.634
7.2 0.546
7.4 0.505
7.5 0.5
lines
10
#100
mappingName
bottomSpline
exit
spline
enter spline points
6
6.5 1.5
6.55 1.5
6.6 1.5
6.7 1.5
6.9 1.5
7. 1.5
lines
10
#100
mappingName
leftSpline
exit
tfi
mappingName
domain3
choose bottom
bottomSpline
choose top
topSpline
choose left
leftSpline
choose right
rightSpline
lines
10 10
#100 100
boundary conditions
1 0 0 1
share
0 0 2 0
exit
check
domain3
#Create Domain 2 ###############
spline
enter spline points
9
6.5 1.5
6.55 1.188
6.6 1.064
6.7 0.9
6.8 0.789
7. 0.634
7.2 0.546
7.4 0.505
7.5 0.5
lines
10
#100
mappingName
rightSplineD2
exit
spline
enter spline points
9
2.5 0.5
2.8 0.5
3.1 0.5
3.4 0.5
4 0.5
5.5 0.5
6 0.5
6.6 0.5
7.5 0.5
lines
10
#100
mappingName
bottomSplineD2
exit
spline
enter spline points
8
2.5 1.5
2.8 1.5
3.1 1.5
3.4 1.5
3.8 1.5
4.2 1.5
5.6 1.5
6.5 1.5
lines
10
#100
mappingName
topSplineD2
exit
spline
enter spline points
9
2.5 1.5
2.5 1.3
2.5 1.2
2.5 1.0
2.5 0.8
2.5 0.7
2.5 0.6
2.5 0.55
2.5 0.5
lines
10
#100
mappingName
leftSplineD2
exit
tfi
mappingName
domain2
choose bottom
bottomSplineD2
choose top
topSplineD2
choose left
leftSplineD2
choose right
rightSplineD2
lines
10 10
#100 100
boundary conditions
0 0 1 0
share
1 2 0 3
exit
check
domain2
#Create Domain 6 ###############
spline
enter spline points
5
13. -0.5
12.875 -0.169
12.75 -0.067
12.625 -0.016
12.5 0.
lines
10
#100
mappingName
bottomSplineRight
exit
spline
enter spline points
7
12.5 0.5
12.5 0.45
12.5 0.4
12.5 0.3
12.5 0.2
12.5 0.1
12.5 0.
lines
10
#100
mappingName
leftSplineRight
exit
spline
enter spline points
9
13.5 -0.5
13.45 -0.188
13.4 -0.0641
13.3 0.1
13.2 0.214
13 0.366
12.8 0.454
12.6 0.495
12.5 0.5
lines
10
#100
mappingName
topSplineRight
exit
spline
enter spline points
6
15.5 -0.5
15.55 -0.5
15.6 -0.5
15.8 -0.5
15.9 -0.5
16 -0.5
lines
10
#100
mappingName
rightSplineRight
exit
tfi
mappingName
domain6
choose bottom
bottomSplineRight
choose top
topSplineRight
choose left
leftSplineRight
choose right
rightSplineRight
lines
10 10
#100 100
boundary conditions
0 1 0 1
exit
check
domain6
#Create Domain 7 ###############
spline
enter spline points
9
13.5 -0.5
13.45 -0.188
13.4 -0.0641
13.3 0.1
13.2 0.214
13 0.366
12.8 0.454
12.6 0.495
12.5 0.5
lines
10
#100
mappingName
leftSplineRightD2
exit
spline
enter spline points
8
13.5 -0.5
14 -0.5
14.5 -0.5
15 -0.5
15.5 -0.5
16 -0.5
16.5 -0.5
17.5 -0.5
line
10
#100
mappingName
bottomSplineRightD2
exit
spline
enter spline points
9
12.5 0.5
13.5 0.5
14 0.5
14.5 0.5
15 0.5
15.5 0.5
16 0.5
16.5 0.5
17.5 0.5
lines
10
#100
mappingName
topSplineRightD2
exit
spline
enter spline points
9
17.5 0.5
17.5 0.3
17.5 0.2
17.5 0.
17.5 -0.1
17.5 -0.3
17.5 -0.4
17.5 -0.45
17.5 -0.5
lines
10
#100
mappingName
rightSplineRightD2
exit
tfi
mappingName
domain7
choose bottom
bottomSplineRightD2
choose top
topSplineRightD2
choose left
leftSplineRightD2
choose right
rightSplineRightD2
lines
10 10
#100 100
boundary conditions
0 0 1 0
exit
check
domain7
#Create Domain 8 ###############
rectangle
mappingName
domain8
set corners
17.5 32.5 -0.5 1.
#17.5 32.5 -0.5 0.5
lines
10 10
#100 100
boundary conditions
0 1 1 1
exit
#Create Domain 9 ###############
rectangle
mappingName
domain9
set corners
12.5 17.5 0.5 1.
#12.5 32.5 0.5 1.
lines
10 10
#100 100
boundary conditions
0 0 0 1
exit
exit
generate an overlapping grid
domain1
domain2
domain3
domain4
domain5
domain6
domain7
domain8
domain9
done
change parameters
prevent hole cutting
domain1
all
domain2
all
domain3
all
domain4
all
domain5
all
#domain6
# all
#domain7
# all
#domain8
# all
#domain9
# all
done
#ghost points
# all
#2 2 2 2
ghost points
domain1
2 2 2 2
ghost points
domain2
2 2 2 2
ghost points
domain3
2 2 2 2
ghost points
domain4
2 2 2 2
ghost points
domain5
2 2 2 2
ghost points
domain6
2 2 2 2
ghost points
domain7
2 2 2 2
ghost points
domain8
2 2 2 2
ghost points
domain9
2 2 2 2
exit
compute overlap
exit
save an overlapping grid
combGridNewMeshOverlap5.hdf
combGridNewMeshOverlap5
exit
________________________________
From: overture-bounce@xxxxxxxxxxxxx<mailto:overture-bounce@xxxxxxxxxxxxx>
<overture-bounce@xxxxxxxxxxxxx<mailto:overture-bounce@xxxxxxxxxxxxx>> on behalf
of Bill Henshaw <wdhenshaw@xxxxxxxxx<mailto:wdhenshaw@xxxxxxxxx>>
Sent: Thursday, June 15, 2017 3:59:07 AM
To: overture@xxxxxxxxxxxxx<mailto:overture@xxxxxxxxxxxxx>
Subject: [overture] Re: Matrix has a row with all zero coefficients
Hi Ohi,
I haven't had a chance to look at what you sent in detail but I suggest that
you
could start by using Overture/tests/tcm3.C to solve Laplace equation with
Dirichlet or Neumann BC's on the grid you created. You could then make
changes to get the BC's that you want.
Regards
Bill
On Wed, Jun 14, 2017 at 6:14 PM, Ohiremen L Dibua
<odibua@xxxxxxxxxxxx<mailto:odibua@xxxxxxxxxxxx>> wrote:
Good evening Bill,
Thank you for being so good about providing this open-source software, and this
forum to ask questions about it. I recently re-meshed a prior set-up by
creating more grids, and implementing what looks to me like proper code for
solving it. I am solving a Laplacian equation. My set up has 9 grids. I get an
error saying
Oges:generateMatrix:ERROR matrix has a row with all zero coefficients
The offending equation is (i1,i2,i3,n,grid)=(10,0,0,0,0).
Below is my code up until the point where an error is thrown, and below that is
the cmd I use to generate the relevant composite grid.
LAPLACIAN CODE
// create and read in a CompositeGrid
CompositeGrid cg;
getFromADataBase(cg,nameOfOGFile);
cg.update();
cg.update(MappedGrid::THEvertex | MappedGrid::THEcenter |
MappedGrid::THEvertexBoundaryNormal);
Ogshow show(nameOfShowFile);
show.saveGeneralComment("Solution to air comb-drive");
show.saveGeneralComment("Test file written on February");
// make a grid function to hold coefficients
Range all;
int stencilSize=int( pow(3,cg.numberOfDimensions())+1.5); //Add 1 for
interpolation equations
realCompositeGridFunction coeff(cg,stencilSize,all,all,all);
realCompositeGridFunction q(cg,all,all,all,2);
coeff.setIsACoefficientMatrix(TRUE,stencilSize);
q=0.;
// create grid functions:
realCompositeGridFunction u, f;
u.link(q,Range(0,0));
f.link(q,Range(1,1));
q.setName("q");
q.setName("u",0);
//q.setName("f",1);
cout << "Number of Grids " << cg.numberOfComponentGrids() << endl;
char buffer[80];
Index I1,I2,I3,I1F,I2F,I3F;
Index Ib1,Ib2,Ib3,Ib1F,Ib2F,Ib3F;
Index Ig1,Ig2,Ig3;
Index Ip1,Ip2,Ip3;
CompositeGridOperators op(cg); // create some differential operators
op.setStencilSize(stencilSize);
coeff.setOperators(op);
cout << "Here " << endl;
coeff=op.laplacianCoefficients(); // get the coefficients for the Laplace
operator
// make shorter names for readability
BCTypes::BCNames dirichlet = BCTypes::dirichlet,
neumann = BCTypes::neumann,
extrapolate = BCTypes::extrapolate,
evenSymmetry = BCTypes::evenSymmetry,
allBoundaries = BCTypes::allBoundaries;
cout << "Here 2 " << endl;
const int ltDomain=1, rtDomain=2, lowDomain=3, upDomain=4;
int grid;
grid=0;
MappedGrid & mg0 = cg[grid];
mg0.update(MappedGrid::THEvertex | MappedGrid::THEcenter |
MappedGrid::THEvertexBoundaryNormal);
mg0.setBoundaryCondition(Start,axis1,ltDomain);
mg0.setBoundaryCondition(End,axis1,rtDomain);
mg0.setBoundaryCondition(Start,axis2,lowDomain);
mg0.setBoundaryCondition(End,axis2,upDomain);
coeff[grid].applyBoundaryConditionCoefficients(0,0,dirichlet,ltDomain);
coeff[grid].applyBoundaryConditionCoefficients(0,0,extrapolate,ltDomain);
coeff[grid].applyBoundaryConditionCoefficients(0,0,dirichlet,lowDomain);
coeff[grid].applyBoundaryConditionCoefficients(0,0,extrapolate,lowDomain);
coeff[grid].applyBoundaryConditionCoefficients(0,0,neumann,upDomain);
coeff[grid].finishBoundaryConditions();
getIndex(mg0.indexRange(),I1,I2,I3);
f[grid](I1,I2,I3)=0;
CMD FOR COMPOSITE GRID
create mappings
#Create Domain 1 ###############
rectangle
mappingName
domain1
set corners
-12.5 2.5 0. 1.5
lines
10 10
boundary conditions
1 0 1 1
exit
#Create Domain 4 ###############
rectangle
mappingName
domain4
set corners
2.5 7.5 0 0.5
#-12.5 7.5 0 0.5
lines
10 10
boundary conditions
0 0 1 0
exit
#Create Domain 5 ###############
rectangle
mappingName
domain5
set corners
7.5 12.5 0 1
lines
10 10
boundary conditions
0 0 1 1
exit
#Create Domain 3 ###############
spline
enter spline points
5
7. 1.5
7.125 1.169
7.25 1.067
7.375 1.016
7.5 1.
lines
10
mappingName
topSpline
exit
spline
enter spline points
7
7.5 1.
7.5 0.9
7.5 0.85
7.5 0.8
7.5 0.7
7.5 0.6
7.5 0.5
lines
10
mappingName
rightSpline
exit
spline
enter spline points
9
6.5 1.5
6.55 1.188
6.6 1.064
6.7 0.9
6.8 0.789
7. 0.634
7.2 0.546
7.4 0.505
7.5 0.5
lines
10
mappingName
bottomSpline
exit
spline
enter spline points
6
6.5 1.5
6.55 1.5
6.6 1.5
6.7 1.5
6.9 1.5
7. 1.5
lines
10
mappingName
leftSpline
exit
tfi
mappingName
domain3
choose bottom
bottomSpline
choose top
topSpline
choose left
leftSpline
choose right
rightSpline
lines
10 10
boundary conditions
1 0 0 1
exit
check
domain3
#Create Domain 2 ###############
spline
enter spline points
9
6.5 1.5
6.55 1.188
6.6 1.064
6.7 0.9
6.8 0.789
7. 0.634
7.2 0.546
7.4 0.505
7.5 0.5
lines
10
mappingName
rightSplineD2
exit
spline
enter spline points
8
2.5 1.5
2.8 1.5
3.1 1.5
3.4 1.5
3.8 1.5
4.2 1.5
5.6 1.5
6.5 1.5
lines
10
mappingName
topSplineD2
exit
spline
enter spline points
9
2.5 0.5
2.8 0.5
3.1 0.5
3.4 0.5
4 0.5
5.5 0.5
6 0.5
6.6 0.5
7.5 0.5
lines
10
mappingName
bottomSplineD2
exit
spline
enter spline points
9
2.5 1.5
2.5 1.3
2.5 1.2
2.5 1.0
2.5 0.8
2.5 0.7
2.5 0.6
2.5 0.55
2.5 0.5
lines
10
mappingName
leftSplineD2
exit
tfi
mappingName
domain2
choose bottom
bottomSplineD2
choose top
topSplineD2
choose left
leftSplineD2
choose right
rightSplineD2
lines
10 10
boundary conditions
0 0 1 0
exit
check
domain2
#Create Domain 6 ###############
spline
enter spline points
5
13. -0.5
12.875 -0.169
12.75 -0.067
12.625 -0.016
12.5 0.
lines
10
mappingName
bottomSplineRight
exit
spline
enter spline points
7
12.5 0.5
12.5 0.45
12.5 0.4
12.5 0.3
12.5 0.2
12.5 0.1
12.5 0.
lines
10
mappingName
leftSplineRight
exit
spline
enter spline points
9
13.5 -0.5
13.45 -0.188
13.4 -0.0641
13.3 0.1
13.2 0.214
13 0.366
12.8 0.454
12.6 0.495
12.5 0.5
lines
10
mappingName
topSplineRight
exit
spline
enter spline points
6
15.5 -0.5
15.55 -0.5
15.6 -0.5
15.8 -0.5
15.9 -0.5
16 -0.5
lines
10
mappingName
rightSplineRight
exit
tfi
mappingName
domain6
choose bottom
bottomSplineRight
choose top
topSplineRight
choose left
leftSplineRight
choose right
rightSplineRight
lines
10 10
boundary conditions
0 1 0 1
exit
check
domain6
#Create Domain 7 ###############
spline
enter spline points
9
13.5 -0.5
13.45 -0.188
13.4 -0.0641
13.3 0.1
13.2 0.214
13 0.366
12.8 0.454
12.6 0.495
12.5 0.5
lines
10
mappingName
leftSplineRightD2
exit
spline
enter spline points
8
13.5 -0.5
14 -0.5
14.5 -0.5
15 -0.5
15.5 -0.5
16 -0.5
16.5 -0.5
17.5 -0.5
line
10
mappingName
bottomSplineRightD2
exit
spline
enter spline points
9
12.5 0.5
13.5 0.5
14 0.5
14.5 0.5
15 0.5
15.5 0.5
16 0.5
16.5 0.5
17.5 0.5
lines
10
mappingName
topSplineRightD2
exit
spline
enter spline points
9
17.5 0.5
17.5 0.3
17.5 0.2
17.5 0.
17.5 -0.1
17.5 -0.3
17.5 -0.4
17.5 -0.45
17.5 -0.5
lines
10
mappingName
rightSplineRightD2
exit
tfi
mappingName
domain7
choose bottom
bottomSplineRightD2
choose top
topSplineRightD2
choose left
leftSplineRightD2
choose right
rightSplineRightD2
lines
10 10
boundary conditions
0 0 1 0
exit
check
domain7
#Create Domain 8 ###############
rectangle
mappingName
domain8
set corners
17.5 32.5 -0.5 1.
#17.5 32.5 -0.5 0.5
lines
10 10
boundary conditions
0 1 1 1
exit
#Create Domain 9 ###############
rectangle
mappingName
domain9
set corners
12.5 17.5 0.5 1.
#12.5 32.5 0.5 1.
lines
10 10
boundary conditions
0 0 0 1
exit
exit
generate an overlapping grid
domain1
domain2
domain3
domain4
domain5
domain6
domain7
domain8
domain9
done
change parameters
prevent hole cutting
#domain1
#all
#domain2
#all
#domain3
#all
#domain4
#all
domain5
all
#domain6
#all
#domain7
#all
#domain8
#all
#domain9
#all
done
ghost points
all
1 1 1 1
exit
compute overlap
exit
save an overlapping grid
combGridNewMeshOverlap5.hdf
combGridNewMeshOverlap5
exit
Thank you for your time,
Ohi
