source: tol/tol_pkg/MulVarFunAprox/DesignC2.Monomial.Spheric.Grid.tol @ 65c7ea3

1941-docker-image-to-testdevel
Last change on this file since 65c7ea3 was 65c7ea3, checked in by Víctor de Buen <vdebuen@…>, 17 years ago

Nuevo parámetro de configuración Real _.radius.margin

git-svn-id: ​https://www.tol-project.org/svn/tolp/trunk@1917 be4846ae-f051-483a-831f-b64e7d08527f
  • Property mode set to 100644
File size: 8.0 KB
Line 
1/* -*- mode: C++ -*- */
2//////////////////////////////////////////////////////////////////////////////
3// FILE   : DesignC2.tol
4// PURPOSE: Defines Class @DesignC2
5//////////////////////////////////////////////////////////////////////////////
6
7
8//////////////////////////////////////////////////////////////////////////////
9Text _.autodoc.member.DesignC2.Monomial.Spheric.Grid =
10"";
11Class @DesignC2.Monomial.Spheric.Grid :
12  @DesignC2.Monomial,
13  @DesignC2.Spheric
14//////////////////////////////////////////////////////////////////////////////
15{
16  @IndexedGrid.Border _.grid;
17//////////////////////////////////////////////////////////////////////////////
18//Configuration parameters
19//////////////////////////////////////////////////////////////////////////////
20  Real _.radius;
21  Real _.grid.size;
22
23//////////////////////////////////////////////////////////////////////////////
24//Auxiliar members
25//////////////////////////////////////////////////////////////////////////////
26  VMatrix _.Y.neighbourhood = Rand(0,0,0,0);
27
28  ////////////////////////////////////////////////////////////////////////////
29  Static @DesignC2.Monomial.Spheric.Grid New(
30    VMatrix X,
31    Real radius.margin,
32    Real radius,
33    Real do.build.border)
34  ////////////////////////////////////////////////////////////////////////////
35  {
36    @DesignC2.Monomial.Spheric.Grid new =
37    [[
38      VMatrix _.X = X;
39      Real    _.n = VColumns(X);
40      Real    _.J = VRows(X);
41      Real    _.radius = radius;
42      Real    _.radius.margin = radius.margin;
43      Real    _.grid.size = radius/(1+radius.margin);
44      @IndexedGrid.Border _.grid = @IndexedGrid.Border::New(
45         X,_.grid.size,True,do.build.border,do.build.border)
46    ]]
47  };
48
49  //////////////////////////////////////////////////////////////////////////////
50  Real build.reference.points(Real void)
51  //////////////////////////////////////////////////////////////////////////////
52  {
53    Text _MID = "[@DesignC2.Monomial.Spheric.Grid::build.reference.points] ";
54    VMatrix _.Y := _.grid::_.Y;
55    VMatrix Y.h = Mat2VMat(_.grid::_.Y.h);
56    VMatrix _.Y.neighbourhood := Y.h | Group("ConcatColumns",
57    For(1,2*_.n, VMatrix(Real i)
58    {
59      VMatrix aux = SubCol(_.grid::_.Y.neighbourhood.index,[[i]]);
60      IfVMat(aux,aux,Y.h)
61    }));
62    Set X.class = Classify(Col2Set(_.grid::point2index(_.X)),Compare);
63    VMatrix _.Yj := Mat2VMat(SetMat([[EvalSet(X.class,Card)]]),1);
64    Real _.X.maxNeighbour := 1+2*_.n;
65    Real _.H := _.grid::_.H
66  };
67
68  //////////////////////////////////////////////////////////////////////////////
69  Real build.stress.points(Real void)
70  //////////////////////////////////////////////////////////////////////////////
71  {
72    Text _MID = "[@DesignC2.Monomial.Spheric.Grid::build.stress.points] ";
73    WriteLn("TRACE "+_MID+"1");
74    Real _.Z.maxNeighbour := If(prior.Z.order>=0,2,0);
75    Matrix Y.neighbour.triplet = Group("ConcatRows", For(1,_.H,Matrix(Real h)
76    {
77      Set neighbours = Unique(VRow2Set(SubRow(_.Y.neighbourhood,[[h]])));
78      Set neighbours.up = Select(neighbours,Real(Real s) { s>h });
79      Matrix s = SetCol(Sort(neighbours.up,Compare));
80      Rand(Rows(s),1,h,h)|s|Rand(Rows(s),1,1,1)
81    }));
82
83    WriteLn("TRACE "+_MID+"2");
84    VMatrix Y.neighbour.1.map = Triplet(Y.neighbour.triplet,_.H,_.H);
85    VMatrix Y.neighbour.2.map = Y.neighbour.1.map*Y.neighbour.1.map;
86    VMatrix Y.neighbour.12.map = Max(Y.neighbour.1.map,Y.neighbour.2.map);
87    Matrix Y.neighbour.12.triplet = VMat2Triplet(Y.neighbour.12.map);
88    VMatrix Z.h.a = Mat2VMat(SubCol(Y.neighbour.12.triplet,[[1]]));
89    VMatrix Z.h.b = Mat2VMat(SubCol(Y.neighbour.12.triplet,[[2]]));
90    VMatrix Y.neighbour.a = SubRow(_.Y,VCol2Set(Z.h.a));
91    VMatrix Y.neighbour.b = SubRow(_.Y,VCol2Set(Z.h.b));
92    VMatrix Z.all = (Y.neighbour.a + Y.neighbour.b)*0.5;
93
94    WriteLn("TRACE "+_MID+"3");
95    VMatrix Z.interior.a =
96    {
97      VMatrix Z.norm2 = RPow(Y.neighbour.a-Z.all,2)*Rand(_.n,1,1,1);
98      VMatrix Z.tao = RPow(Z.norm2,0.5)/_.radius;
99      LE(Z.tao,Rand(VRows(Z.tao),1,1,1))
100    };
101    VMatrix Z.interior.b =
102    {
103      VMatrix Z.norm2 = RPow(Y.neighbour.b-Z.all,2)*Rand(_.n,1,1,1);
104      VMatrix Z.tao = RPow(Z.norm2,0.5)/_.radius;
105      LE(Z.tao,Rand(VRows(Z.tao),1,1,1))
106    };
107    VMatrix match.interior = Z.interior.a $* Z.interior.b;
108
109    Set sel.interior = MatQuery::SelectMatch(match.interior);
110    VMatrix Z.interior = SubRow(Z.all,sel.interior);
111    VMatrix Z.h.interior = SubRow(Z.h.a | Z.h.b,sel.interior);
112    Set criterium = Range(1,_.n,1);
113    Set sel.unique = Col2Set(Sort(VMat2Mat(Z.interior),criterium,True));
114    WriteLn("TRACE "+_MID+"4");
115    VMatrix _.Z := SubRow(Z.interior, sel.unique);
116    VMatrix _.Z.h := SubRow(Z.h.interior, sel.unique);
117    Real _.Tz := VRows(_.Z);
118    VMatrix Z.Y.neighbour =
119    {
120      Matrix Z.t = DifEq(1/(1-B),Rand(_.Tz,1,1,1));
121      Matrix Z.th = (Z.t << Z.t) |
122       VMat2Mat((SubCol(_.Z.h,[[1]])<<SubCol(_.Z.h,[[2]])));
123      Matrix Z.Y.idx = PivotByRows(Z.th,Sort(Z.th,[[2]]));
124      VMatrix Z.Y.dif =
125       SubRow(_.Z,Col2Set(SubCol(Z.Y.idx,[[1]]))) -
126       SubRow(_.Y,Col2Set(SubCol(Z.Y.idx,[[2]])));
127        VMatrix Z.Y.norm2 = RPow(Z.Y.dif,2)*Rand(_.n,1,1,1);
128        VMatrix Z.Y.tao = RPow(Z.Y.norm2,0.5)/_.radius;
129        VMatrix Z.Y.match = LE(Z.Y.tao,Rand(VRows(Z.Y.tao),1,1,1));
130      SubRow(Mat2VMat(Z.Y.idx),MatQuery::SelectMatch(Z.Y.match))
131    };
132    VMatrix aux = MatQuery::SortAndCount(SubCol(Z.Y.neighbour,[[2]]));
133    VMatrix _.Yt := Triplet(VMat2Mat(SubCol(aux,[[1]])) | 
134                            Rand(VRows(aux),1,1,1) |
135                            VMat2Mat(SubCol(aux,[[2]])),_.H,1);
136    Real _.Mz := Case(
137      prior.Z.order< 0, 0,
138      prior.Z.order==0, _.Tz,
139      prior.Z.order==1, _.Tz*(1+_.n),
140      prior.Z.order==2, _.Tz*(1+_.n*3/2+(1/2)*_.n^2));
141    _.Tz
142  };
143
144  //////////////////////////////////////////////////////////////////////////////
145  Real build.border.points(Real void)
146  //////////////////////////////////////////////////////////////////////////////
147  {
148    If(!Rows(_.grid::_.border.assign),{
149      VMatrix _.Yb := Zeros(_.H,1);
150      Real _.Tb := 0
151    },{
152    Text _MID = "[@DesignC2.Monomial.Spheric.Grid::build.border.points] ";
153    VMatrix   B.t  = Mat2VMat(SubCol(_.grid::_.border.assign,[[1]]));
154    VMatrix _.B.h := Mat2VMat(SubCol(_.grid::_.border.assign,[[2]]));
155    VMatrix _.B := SubRow(_.grid::_.border::_.Y, VCol2Set(B.t));
156    VMatrix aux = MatQuery::SortAndCount(_.B.h);
157    VMatrix _.Yb := Triplet(VMat2Mat(SubCol(aux,[[1]])) | 
158                            Rand(VRows(aux),1,1,1) |
159                            VMat2Mat(SubCol(aux,[[2]])),_.H,1);
160
161    Real _.Tb := VRows(_.B);
162    Real _.Mb := Case(
163      prior.B.order< 0, 0,
164      prior.B.order==0, _.Tb,
165      prior.B.order==1, _.Tb*(1+_.n),
166      prior.B.order==2, _.Tb*(1+_.n*3/2+(1/2)*_.n^2));
167    _.Tb
168    })
169  };
170
171  //////////////////////////////////////////////////////////////////////////////
172  VMatrix get.radius(VMatrix Y.idx)
173  //////////////////////////////////////////////////////////////////////////////
174  {
175    Constant(VRows(Y.idx),1,_.radius)
176  };
177
178  //////////////////////////////////////////////////////////////////////////////
179  VMatrix closest.reference.index(VMatrix X)
180  //////////////////////////////////////////////////////////////////////////////
181  {
182    WriteLn("TRACE [closest.reference.index] 1");
183    SubRow(_.Y.neighbourhood, Col2Set(_.grid::point2index(X)))
184  };
185
186  ////////////////////////////////////////////////////////////////////////////
187    VMatrix match.evaluable(VMatrix X)
188  ////////////////////////////////////////////////////////////////////////////
189  {
190    Real m = VRows(X);
191    Mat2VMat(_.grid::point2index(X))
192  };
193
194  ///////////////////////////////////////////////////////////////////////////
195    VMatrix generate.internal.points(Real m)
196  ///////////////////////////////////////////////////////////////////////////
197  {
198    VMatrix h = Round(Rand(m,1,0.5,_.H+0.499999999));
199    VMatrix y = SubRow(_.Y,VCol2Set(h));
200    VMatrix v = Gaussian(m,_.n,0,1);
201    VMatrix u = v*(VMatDat(MtMSqr(v),1,1)^-0.5);
202    y+u*(Rand(0,1)*_.radius)
203  }
204
205
206};
Note: See TracBrowser for help on using the repository browser.