Changeset b050cc43 in tol for tol_pkg/MulVarFunAprox


Ignore:
Timestamp:
Dec 4, 2009, 4:44:16 PM (17 years ago)
Author:
Víctor de Buen <vdebuen@…>
Branches:
1941-docker-image-to-test, devel, master
Children:
4bd88af
Parents:
fabac4a
Message:

New class hierarchy

git-svn-id: ​https://www.tol-project.org/svn/tolp/trunk@1762 be4846ae-f051-483a-831f-b64e7d08527f

Location:
tol_pkg/MulVarFunAprox
Files:
4 added
5 edited

Legend:

Unmodified
Added
Removed
  • tol_pkg/MulVarFunAprox/DesignC2.Monomial.Spheric.Grid.tol

    rfabac4a rb050cc43  
    99Text _.autodoc.member.DesignC2.Monomial.Spheric.Grid =
    1010"";
    11 Class @DesignC2.Monomial.Spheric.Grid : @DesignC2
     11Class @DesignC2.Monomial.Spheric.Grid :
     12  @DesignC2.Monomial,
     13  @DesignC2.Spheric
    1214//////////////////////////////////////////////////////////////////////////////
    1315{
    … …  
    1719 
    1820  Real force.reference.near.nodal = True;
    19   Real radius.margin  = 0.05;
    20   Real max.deg = 6;
    21   Real grid.size = ?;
     21  Real _.radius = ?;
     22  Real _.grid.size = ?;
    2223
    2324//////////////////////////////////////////////////////////////////////////////
    … …  
    2728  Set _.X.col.info = Copy(Empty);
    2829  Set _.range = Copy(Empty);
    29   Real _.radius = ?;
    3030  VMatrix _.grid = Rand(0,0,0,0);
    3131  Real _.G = ?;
    … …  
    3737
    3838  ////////////////////////////////////////////////////////////////////////////
    39   Static @DesignC2.Monomial.Spheric.Grid New(VMatrix X)
     39  Static @DesignC2.Monomial.Spheric.Grid New(VMatrix X, Real radius.margin)
    4040  ////////////////////////////////////////////////////////////////////////////
    4141  {
    … …  
    4444      VMatrix _.X = X;
    4545      Real    _.n = VColumns(X);
    46       Real    _.J = VRows(X)
     46      Real    _.J = VRows(X);
     47      Real _.radius.margin = radius.margin
    4748    ]]
    4849  };
    … …  
    6162      VMatrix Xi = SubCol(X,[[i]]);
    6263      VMatrix Yi.min = Rand(m,1,_.range[i][1],_.range[i][1]);
    63       VMatrix Ri = Round((Xi-Yi.min)/grid.size);
     64      VMatrix Ri = Round((Xi-Yi.min)/_.grid.size);
    6465      VMatrix minRi = Rand(m,1,0,0);
    6566      VMatrix maxRi = Rand(m,1,cRi-1,cRi-1);
    … …  
    99100    ]]});
    100101   
    101     Real If(Or(IsUnknown(grid.size),grid.size<=0),
     102    Real If(Or(IsUnknown(_.grid.size),_.grid.size<=0),
    102103    {
    103104      Real length.max = SetMax(Traspose(Extract(_.X.col.info,4))[1]);
    104105      Real density = _.J/SetProd(Traspose(Extract(_.X.col.info,4))[1]);
    105       grid.size := 0.5*(_.K0/density)^(1/_.n)
    106     });
    107    
    108     Real _.radius := (1+radius.margin)*grid.size/Sqrt(2);
     106      _.grid.size := 0.5*(_.K0/density)^(1/_.n)
     107    });
     108   
     109    Real _.radius := (1+_.radius.margin)*_.grid.size/Sqrt(2);
    109110   
    110111    WriteLn("TRACE "+_MID+"2");
    … …  
    112113    {
    113114      Set info = _.X.col.info[i];
    114       Real points = 1+Ceil(info::length/grid.size);
    115       Real exceed = (points-1)*grid.size-info::length;
     115      Real points = 1+Ceil(info::length/_.grid.size);
     116      Real exceed = (points-1)*_.grid.size-info::length;
    116117      Matrix rnk = Tra(DifEq(1/(1-B),Rand(points,1,1,1))-Rand(points,1,1,1));   
    117       Set MatSet((rnk*grid.size)+info::min-exceed/2)[1]
     118      Set MatSet((rnk*_.grid.size)+info::min-exceed/2)[1]
    118119    });
    119120    WriteLn("TRACE "+_MID+"3");
    … …  
    165166    VMatrix Y.closest.grid.index = closest.grid.index(_.Y);
    166167   
    167     VMatrix Y.neighbour.pos = Min(Y.max, _.Y + Rand(_.H,_.n,grid.size,grid.size));
    168     VMatrix Y.neighbour.neg = Max(Y.min, _.Y - Rand(_.H,_.n,grid.size,grid.size));
     168    VMatrix Y.neighbour.pos = Min(Y.max, _.Y + Rand(_.H,_.n,_.grid.size,_.grid.size));
     169    VMatrix Y.neighbour.neg = Max(Y.min, _.Y - Rand(_.H,_.n,_.grid.size,_.grid.size));
    169170   
    170171    VMatrix _.Y.h := DifEq(1/(1-B),Rand(_.H,1,1,1));
    … …  
    249250    VMatrix _.Z.h := SubRow(Z.h.a | Z.h.b,MatQuery::SelectMatch(Z.interior));
    250251   
    251     Real _.T := VRows(_.Z);
     252    Real _.Tz := VRows(_.Z);
    252253    Real _.Mz := Case(
    253254      prior.Z.order< 0, 0,
    254       prior.Z.order==0, _.T,
    255       prior.Z.order==1, _.T*(1+_.n),
    256       prior.Z.order==2, _.T*(1+_.n*3/2+(1/2)*_.n^2));
    257     _.T
    258   };
    259 
    260   //////////////////////////////////////////////////////////////////////////////
    261     Real build.base(Real void)
    262   //////////////////////////////////////////////////////////////////////////////
    263   {
    264     Text _MID = "[@DesignC2.Monomial.Spheric.Grid::build.base] ";
    265     Real max.K = Min(_.K0, Floor(_.M/(_.H*min.contrast.surface)));
    266    
    267     WriteLn("TRACE "+_MID+"1.1 K0="<<_.K0);
    268     WriteLn("TRACE "+_MID+"1.2 M="<<_.M);
    269     WriteLn("TRACE "+_MID+"1.3 H="<<_.H);
    270     WriteLn("TRACE "+_MID+"1.4 SC="<<min.contrast.surface);
    271     WriteLn("TRACE "+_MID+"1.5 max.K="<<max.K);
    272     Real deg =
    273     {
    274       Real deg = 0;
    275       Real K = 1;
    276       Set pair.deg.K = [[ [[Copy(deg),Copy(K),1]] ]];
    277       While(And(deg<=max.deg, K<=max.K),
    278       {
    279         Real deg := deg+1;
    280         Real k = Comb(_.n+deg-1,deg);
    281         Real K := K+k;
    282         WriteLn("TRACE "+_MID+"2."<<deg+" K="<<K);
    283         Set Append(pair.deg.K, [[ [[Copy(deg),Copy(K),Copy(k)]] ]])
    284       });
    285       If(pair.deg.K[deg+1][2]<=max.K,deg,deg-1)
    286     };
    287 
    288     WriteLn("TRACE "+_MID+"3");
    289     FunRn2R::@BaseC2 base.feasible = FunRn2R::@BaseC2.Monomial::New(_.n,deg);
    290     Real _.K := base.feasible::_.K;
    291     Set Append(_.base, [[base.feasible]]);
    292     Real _.base.sel := Card(_.base);
    293     _.K
     255      prior.Z.order==0, _.Tz,
     256      prior.Z.order==1, _.Tz*(1+_.n),
     257      prior.Z.order==2, _.Tz*(1+_.n*3/2+(1/2)*_.n^2));
     258    _.Tz
    294259  };
    295260
    … …  
    305270    Real _.K0 := base.max::_.K;
    306271    _build.model.design(void)
    307   };
    308 
    309   ///////////////////////////////////////////////////////////////////////////
    310   Set get.X.weight(Set X.dif)
    311   ///////////////////////////////////////////////////////////////////////////
    312   {
    313     Real m = VRows(X.dif[1]);
    314     Set X.tao = For(1,_.X.maxNeighbour,VMatrix(Real i)
    315     {
    316       VMatrix X.norm2 = RPow(X.dif[i],2)*Rand(_.n,1,1,1);
    317       RPow(X.norm2,0.5)/_.radius
    318     });
    319     Set X.tao2 = For(1,_.X.maxNeighbour,VMatrix(Real i)
    320     {
    321       RPow(X.tao[i],2)
    322     });
    323     Set X.match = For(1,_.X.maxNeighbour,VMatrix(Real i)
    324     {
    325       LE(X.tao[i],Rand(m,1,1,1))
    326     });
    327     Set X.lambda = For(1,_.X.maxNeighbour,VMatrix(Real i)
    328     {
    329       RPow(X.match[i] $* (Rand(m,1,1,1)-X.tao2[i]),4)
    330     });
    331     VMatrix X.lambda.sum = SetSum(X.lambda);
    332    
    333     Set X.weight = For(1,_.X.maxNeighbour,VMatrix(Real i)
    334     {
    335       X.lambda[i] $/ X.lambda.sum
    336     })
    337272  };
    338273
  • tol_pkg/MulVarFunAprox/DesignC2.tol

    rfabac4a rb050cc43  
    2828  Real _.X.maxNeighbour = ?;
    2929  Real _.Z.maxNeighbour = ?;
    30   VMatrix _.Y = Rand(0,0,0,0);
    31   Real _.H = ?;
    32   VMatrix _.Z = Rand(0,0,0,0);
    33   VMatrix _.Z.h = Rand(0,0,0,0);
    34   Real _.T = ?;
    35   Real _.Mz = ?;
    36   Real _.M = ?;
    37   Real _.K = ?;
    38   Real _.N = ?;
     30  VMatrix _.Y = Rand(0,2,0,0);
     31  VMatrix _.Z = Rand(0,2,0,0);
     32  VMatrix _.Z.h = Rand(0,1,0,0);
     33  VMatrix _.B = Rand(0,2,0,0);
     34  VMatrix _.B.h = Rand(0,0,0,0);
     35  Real _.K = 0;
     36  Real _.H = 0;
     37  Real _.Tz = 0;
     38  Real _.Mz = 0;
     39  Real _.Tb = 0;
     40  Real _.Mb = 0;
     41  Real _.M = 0;
     42  Real _.N = 0;
    3943  Set _.base = Copy(Empty);
    4044  Real _.base.sel = ?;
    4145  VMatrix _.QX = Rand(0,0,0,0);
     46  VMatrix _.QB = Rand(0,0,0,0);
    4247  VMatrix _.PZ = Rand(0,0,0,0);
    4348  VMatrix _.Q = Rand(0,0,0,0);
    4449  VMatrix _.LQ = Rand(0,0,0,0);
    45 
     50  VMatrix _.Q.nonNullCellsByColumns = Rand(0,0,0,0);
    4651  //////////////////////////////////////////////////////////////////////////////
    4752  NameBlock get.base(Real sel)
    … …  
    5560  //Virtual methods
    5661  //////////////////////////////////////////////////////////////////////////////
     62  Real build.base(Real void);
    5763  Real build.reference.points(Real void);
    58   Real build.stress.points(Real void);
    59   Real build.base(Real void);
     64  Real build.stress.points(Real void) { False };
     65  Real build.border.points(Real void) { False };
     66  VMatrix closest.reference.index(VMatrix X);
    6067  Set get.X.weight(Set X.dif);
    6168  VMatrix match.evaluable(VMatrix X);
    62   VMatrix closest.reference.index(VMatrix X);
    6369
    6470  //////////////////////////////////////////////////////////////////////////////
    … …  
    7076    Real build.reference.points(0);
    7177    Real build.stress.points(0);
    72     Real _.M := _.J+_.Mz;
     78    Real build.border.points(0);
     79    Real _.M := _.J+_.Mb+_.Mz;
    7380    Real build.base(0);
    7481    Real _.N := _.H*_.K;
    … …  
    8491
    8592  //////////////////////////////////////////////////////////////////////////////
    86     VMatrix build.obs(VMatrix X)
    87   //////////////////////////////////////////////////////////////////////////////
    88   {
     93    VMatrix build.obs(VMatrix X, VMatrix X.h)
     94  //////////////////////////////////////////////////////////////////////////////
     95  {
     96    Text _MID = "[@DesignC2::build.obs] ";
    8997    Real m = VRows(X);
    90     VMatrix X.h=closest.reference.index(X);
     98    WriteLn("TRACE "+_MID+"1 m="<<m);
     99    If(!m, Rand(0,_.N,0,0), {
     100    WriteLn("TRACE "+_MID+"2");
     101    WriteLn("TRACE "+_MID+"2");
    91102    VMatrix X.j=DifEq(1/(1-B),Rand(m,1,1,1));
     103    WriteLn("TRACE "+_MID+"2");
    92104    Set X.dif = For(1,_.X.maxNeighbour,VMatrix(Real i)
    93105    {
    … …  
    95107      X - SubRow(_.Y,X.h.i)
    96108    });
     109    WriteLn("TRACE "+_MID+"2");
    97110    Set X.weight = get.X.weight(X.dif);
     111    WriteLn("TRACE "+_MID+"2");
    98112    FunRn2R::@BaseC2 b = get.base(0);
     113    WriteLn("TRACE "+_MID+"2");
    99114    Matrix QX.triplet = Group("ConcatRows", For(1,_.X.maxNeighbour,Matrix(Real i)
    100115    {
    … …  
    114129    }));
    115130    Pack(Convert(Triplet(QX.triplet, m, _.N),"Cholmod.R.Sparse"))
    116   };
    117   //////////////////////////////////////////////////////////////////////////////
    118     VMatrix build.prior(Real void)
    119   //////////////////////////////////////////////////////////////////////////////
    120   {
    121     VMatrix Z.t=DifEq(1/(1-B),Rand(_.T,1,1,1));
     131  })};
     132  //////////////////////////////////////////////////////////////////////////////
     133    VMatrix build.border.prior(Real void)
     134  //////////////////////////////////////////////////////////////////////////////
     135  {
     136    If(!_.Tb, Rand(0,_.N,0,0), {
     137    VMatrix B.t=DifEq(1/(1-B),Rand(_.Tb,1,1,1));
     138    VMatrix Y = SubRow(_.Y, VCol2Set(_.B.h));
     139    VMatrix B.dif = _.B - Y;
     140    FunRn2R::@BaseC2 b = get.base(0);   
     141    Matrix PB.triplet =
     142    {
     143      VMatrix B.hK = (_.B.h-1)*_.K;
     144      Real b::upgrade.cache(B.dif,prior.Z.order);
     145      Real lastRow = 0;
     146      Group("ConcatRows", For(1,_.K,Matrix(Real k)
     147      {
     148        Real cumRows = lastRow;
     149        VMatrix pk = b::function.v(k,B.dif);
     150        VMatrix col.idx = (B.hK+k);
     151        {
     152          VMatrix row.idx = B.t+cumRows;
     153          Real cumRows := cumRows + _.Tb;
     154          VMat2Mat(row.idx | col.idx | pk)
     155        }  << If(prior.Z.order<1,Rand(0,3,0,0),
     156        Group("ConcatRows",For(1,_.n,Matrix(Real d)
     157        {
     158          VMatrix pkd = b::partial.derivative.1.v(k,B.dif,d);
     159          VMatrix row.idx = B.t+cumRows;
     160          Real cumRows := cumRows + _.Tb;
     161          VMat2Mat(row.idx | col.idx | pkd)
     162        })))<< If(prior.Z.order<2,Rand(0,3,0,0),
     163        Group("ConcatRows",For(1,_.n,Matrix(Real d1)
     164        {
     165          Group("ConcatRows",For(1,d1,Matrix(Real d2)
     166          {
     167            VMatrix pkd2 = b::partial.derivative.2.v(k,B.dif,d1,d2);
     168            VMatrix row.idx = B.t+cumRows;
     169            Real cumRows := cumRows + _.Tb;
     170            VMat2Mat(row.idx | col.idx | pkd2)
     171          }))
     172        })))
     173      }))
     174    };
     175    Convert(Triplet(PB.triplet, _.Mb, _.N),"Cholmod.R.Sparse")
     176  })};
     177  //////////////////////////////////////////////////////////////////////////////
     178    VMatrix build.stress.prior(Real void)
     179  //////////////////////////////////////////////////////////////////////////////
     180  {
     181    If(!_.Tz, Rand(0,_.N,0,0), {
     182    VMatrix Z.t=DifEq(1/(1-B),Rand(_.Tz,1,1,1));
    122183   
    123184    Set Z.dif = For(1,_.Z.maxNeighbour,VMatrix(Real i)
    … …  
    140201        {
    141202          VMatrix row.idx = Z.t.i+cumRows;
    142           Real cumRows := cumRows + _.T;
     203          Real cumRows := cumRows + _.Tz;
    143204          VMat2Mat(row.idx | col.idx | pk)
    144205        }  << If(prior.Z.order<1,Rand(0,3,0,0),
    … …  
    148209          VMatrix pkd = b::partial.derivative.1.v(k,Zdif.i,d)*sign;
    149210          VMatrix row.idx = Z.t.i+cumRows;
    150           Real cumRows := cumRows + _.T;
     211          Real cumRows := cumRows + _.Tz;
    151212          VMat2Mat(row.idx | col.idx | pkd)
    152213        })))<< If(prior.Z.order<2,Rand(0,3,0,0),
    … …  
    158219            VMatrix pkd2 = b::partial.derivative.2.v(k,Zdif.i,d1,d2)*sign;
    159220            VMatrix row.idx = Z.t.i+cumRows;
    160             Real cumRows := cumRows + _.T;
     221            Real cumRows := cumRows + _.Tz;
    161222            VMat2Mat(row.idx | col.idx | pkd2)
    162223          }))
    … …  
    165226    }));
    166227    Convert(Triplet(PZ.triplet, _.Mz, _.N),"Cholmod.R.Sparse")
    167   };
     228  })};
    168229  //////////////////////////////////////////////////////////////////////////////
    169230    Real build.decomposition(Real void)
    170231  //////////////////////////////////////////////////////////////////////////////
    171232  {
    172     VMatrix _.QX := build.obs(_.X);
    173     VMatrix _.PZ := build.prior(0);
    174     VMatrix Q.all = _.QX << _.PZ;
    175     Set Q.active = MatQuery::SelectMatch(
    176       NE(Tra(Rand(1,_.M,1,1)*Abs(Q.all)),Rand(_.N,1,0,0)));
    177     VMatrix _.Q := SubCol(Q.all,Q.active);
     233    Text _MID = "[@DesignC2::build.decomposition] ";
     234    WriteLn("TRACE "+_MID+"1");
     235    VMatrix X.h=closest.reference.index(_.X);
     236    VMatrix _.QX := build.obs(_.X, X.h);
     237    WriteLn("TRACE "+_MID+"2");
     238    VMatrix _.QB := build.border.prior(0);
     239    WriteLn("TRACE "+_MID+"3");
     240    VMatrix _.PZ := build.stress.prior(0);
     241    WriteLn("TRACE "+_MID+"4");
     242    VMatrix _.Q := _.QX << _.QB << _.PZ;
     243    WriteLn("TRACE "+_MID+"7");
     244    VMatrix _.Q.nonNullCellsByColumns :=
     245    {
     246      VMatrix aux = Tra(Rand(1,_.M,1,1)*Or(_.Q,Zeros(_.M,_.N)));
     247      Group("ConcatColumns",For(1,_.K,VMatrix(Real k)
     248      {
     249        SubRow(aux,Range(k,_.N,_.K))
     250      }))
     251    };
     252    WriteLn("TRACE "+_MID+"8");
    178253    VMatrix _.LQ := CholeskiFactor(_.Q,"XtX");
     254    WriteLn("TRACE "+_MID+"9");
    179255    True
    180256  }
  • tol_pkg/MulVarFunAprox/ModelC2.tol

    rfabac4a rb050cc43  
    1717//Auxiliar members
    1818//////////////////////////////////////////////////////////////////////////////
    19   VMatrix _.Fz = Rand(0,0,0,0);
    20   VMatrix _.F = Rand(0,0,0,0);
    21   VMatrix _.alpha = Rand(0,0,0,0);
    22   VMatrix _.F.frc = Rand(0,0,0,0);
    23   VMatrix _.F.cmp = Rand(0,0,0,0);
    24   VMatrix _.E = Rand(0,0,0,0);
     19  VMatrix _.F = Rand(0,1,0,0);
     20  VMatrix _.Fz = Rand(0,1,0,0);
     21  VMatrix _.Fb = Rand(0,1,0,0);
     22  VMatrix _.alpha = Rand(0,1,0,0);
     23  VMatrix _.F.frc = Rand(0,1,0,0);
     24  VMatrix _.F.cmp = Rand(0,2,0,0);
     25  VMatrix _.E = Rand(0,1,0,0);
    2526  Real _.E.avr = ?;
    2627  Real _.E.ste = ?;
    … …  
    4445  {
    4546    VMatrix _.Fz := Rand(_.design::_.Mz,1,0,0);
    46     VMatrix _.F := _.Fx<<_.Fz;
     47    VMatrix _.Fb := Rand(_.design::_.Mb,1,0,0);
     48    VMatrix _.F := _.Fx << _.Fb << _.Fz;
    4749    VMatrix _.alpha := CholeskiSolve(_.design::_.LQ,Tra(Tra(_.F)*_.design::_.Q),"PtLLtP");
    4850    VMatrix _.F.frc := SubRow(_.design::_.Q*_.alpha,Range(1,_.design::_.J,1));
    … …  
    6365    VMatrix X.1 = SubRow(X,X.select);
    6466    Real m.1 = VRows(X.1);
    65     VMatrix QX.1 = _.design::build.obs(X.1);
     67    VMatrix X.1.h = _.design::closest.reference.index(X.1);
     68    VMatrix QX.1 = _.design::build.obs(X.1, X.1.h);
    6669    VMatrix F.1 = QX.1*_.alpha;
    6770    Matrix triplet.1 = SetCol(X.select)|Rand(m.1,1,1,1)|VMat2Mat(F.1);
  • tol_pkg/MulVarFunAprox/MulVarFunAprox.tol

    rfabac4a rb050cc43  
    1616   "vdebuen@tol-project.org"   ]];
    1717
     18#Embed "IndexedGrid.tol";
    1819#Embed "DesignC2.tol";
    1920#Embed "DesignC2.Monomial.tol";
    2021#Embed "DesignC2.Spheric.tol";
    2122#Embed "DesignC2.Monomial.Spheric.Grid.tol";
     23#Embed "DesignC2.Monomial.Spheric.Disconnected.tol";
    2224#Embed "ModelC2.tol"
    2325
  • tol_pkg/MulVarFunAprox/test/test_0001/test.tol

    rfabac4a rb050cc43  
    1414
    1515//Random seed setting
    16 Real PutRandomSeed(0);
     16Real PutRandomSeed(1096004697);
    1717Real rndSeed = GetRandomSeed(0);
    1818WriteLn("Current Random Seed = "<<rndSeed);
    … …  
    3333  }))
    3434};
    35 
    36 MulVarFunAprox::@DesignC2 design =
     35VMatrix X = generate.nodal.points(15000,2);
     36/* * /
     37MulVarFunAprox::@DesignC2 design = {
    3738  MulVarFunAprox::@DesignC2.Monomial.Spheric.Grid::New(
    38     generate.nodal.points(15000,2));
    39 
    40 Real design::force.reference.near.nodal := True;
    41 Real design::radius.margin := 0.05;
     39    X, radius.margin=0.05) };
     40/* */
     41MulVarFunAprox::@DesignC2 design = {
     42  MulVarFunAprox::@DesignC2.Monomial.Spheric.Disconnected::New(
     43    VMatrix X,
     44    Real radius.margin = Sqrt(2)-1,
     45    Real radius = .5,
     46    Real do.build.border = False)};
     47/* */
     48//Real design::force.reference.near.nodal := True;
     49Real design::_.radius.margin := 0.05;
    4250Real design::max.deg := 6;
    4351
    … …  
    8795VMatrix X2 = { design::generate.internal.points(
    8896  Real m = design::_.J/2,
    89   Real design::force.reference.near.nodal)};
     97  True)};
    9098
    9199Real b::upgrade.cache(X2,0);
Note: See TracChangeset for help on using the changeset viewer.