source: tol/tol_pkg/MulVarFunAprox/DesignC2.Monomial.Spheric.Grid.tol @ 9e84556

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

Generalización para entornos de radio variable

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