atlas.h File Reference

Detailed Description

Definition of a set of local charts on a manifold.

See Also
Tatlas,atlas.c

Definition in file atlas.h.

Data Structures

struct  TAtlasStatistics
 Data about the atlas construction. More...
 
struct  Tbifurcation
 Representation of a bifurcation. More...
 
struct  TAtlasHeapElement
 Structure used for the elements of the heap of charts. More...
 
struct  Tatlas
 A atlas on a manifold. More...
 

Macros

#define ATLAS_VERBOSE   0
 Vebosity of the atlas operations. More...
 
#define PLOT_ATLAS_IN_COLORS   1
 Activate/deactivate the colors in the atlas plot. More...
 
#define PROJECT_WITH_PLANE   1
 Select between the two types of projection when interpolating between configurations. More...
 
#define INIT_NUM_CHARTS   100
 Initial number of charts of an atlas. More...
 
#define INIT_NUM_BIFURCATIONS   100
 Initial number of bifurcations of an atlas. More...
 
#define USE_ATLAS_TREE   0
 Whether to use a binary tree to search for neighbouring charts. More...
 
#define SAMPLING_RADIUS_REDUCTION_FACTOR   0.9
 Ratio at which the sampling radius is reduced. More...
 
#define SIMPLE_BORDER   1
 If 1, the border of the atlas are roughly approximated. More...
 
#define INIT_NUM_CHARTS_IN_ATLAS_HEAP   100
 Initial number of elements in the heap of local charts. More...
 

Functions

void InitAtlasHeapElement (unsigned int mID, double c, double beta, TAtlasHeapElement *he)
 Constructor of TAtlasHeapElement. More...
 
void CopyAtlasHeapElement (void *he1, void *he2)
 Copy constructor. More...
 
unsigned int GetAtlasHeapElementID (TAtlasHeapElement *he)
 Gets the chart identifier. More...
 
double GetAtlasHeapElementCost (TAtlasHeapElement *he)
 Gets the cost. More...
 
boolean LessThanAtlasHeapElement (void *he1, void *he2, void *userData)
 Comparison between atlas heap elements. More...
 
void PenalizeAtlasHeapElement (TAtlasHeapElement *he)
 Penalized the cost stored in a atlas heap element. More...
 
void DeleteAtlasHeapElement (void *he)
 Atlas heap element destructor. More...
 
void InitAtlas (Tparameters *pr, boolean parallel, boolean simpleChart, unsigned int k, double e, double ce, double r, TAtlasBase *w, Tatlas *a)
 Constructor. More...
 
unsigned int AddChart2Atlas (Tparameters *pr, double *ps, unsigned int parentID, boolean *singularity, Tatlas *a)
 Defines a new chart and adds it to the atlas. More...
 
unsigned int AddTrustedChart2Atlas (Tparameters *pr, double *ps, unsigned int parentID, boolean *singularity, Tatlas *a)
 Defines a new chart and adds it to the atlas. More...
 
boolean InitAtlasFromPoint (Tparameters *pr, boolean parallel, boolean simpleChart, double *p, TAtlasBase *w, Tatlas *a)
 Initializes an atlas from a given point. More...
 
void BuildAtlasFromPoint (Tparameters *pr, double *p, boolean simpleChart, TAtlasBase *w, Tatlas *a)
 Defines an atlas from a given point. More...
 
boolean MinimizeOnAtlas (Tparameters *pr, char *fname, double *p, TAtlasBase *w, double(*costF)(Tparameters *, boolean, double *, void *), void *costData, Tatlas *a)
 Gradient minimization on an manifold. More...
 
double GetAtlasRadius (Tatlas *a)
 Radius used in the charts. More...
 
double GetAtlasError (Tatlas *a)
 Maximum error for all the charts. More...
 
double GetAtlasErrorCurv (Tatlas *a)
 Maximum curvature error for all the charts. More...
 
TAtlasBaseGetAtlasWorld (Tatlas *a)
 Gets the world structure on which the atlas is defined. More...
 
unsigned int GetAtlasAmbientDim (Tatlas *a)
 Ambient dimensionality. More...
 
unsigned int GetAtlasManifoldDim (Tatlas *a)
 Manifold dimensionality. More...
 
boolean Newton2ManifoldPlane (Tparameters *pr, double *point, double *vector, double *pInit, double *p, Tatlas *a)
 Finds a point in the intersection of the manifold and a plane. More...
 
boolean ExtendAtlasFromPoint (Tparameters *pr, unsigned int parentID, double *t, TAtlasStatistics *st, Tatlas *a)
 Tries to expand the atlas from a chart from a given point. More...
 
boolean ExtendAtlasTowardPoint (Tparameters *pr, unsigned int parentID, double *t, boolean collisionStops, TAtlasStatistics *st, Tatlas *a)
 Tries to expand the atlas from a chart toward a given point. More...
 
boolean AtlasAStar (Tparameters *pr, double *p, double *time, double *pl, unsigned int *ns, double ***path, Tatlas *a)
 Expands the atlas to reach a given point. More...
 
boolean AtlasGBF (Tparameters *pr, double *p, double *time, double *pl, unsigned int *ns, double ***path, Tatlas *a)
 Expands the atlas to reach a given point. More...
 
void AtlasTree (Tparameters *pr, unsigned int nNodes, Tatlas *a)
 Defines a tree of charts on a manifold. More...
 
unsigned int GetAtlasNumCharts (Tatlas *a)
 Number of charts in the atlas. More...
 
TchartGetAtlasChart (unsigned int id, Tatlas *a)
 Gets one of the charts of the chart. More...
 
boolean RandomPointInAtlas (Tparameters *pr, double scale, double *w, unsigned int *nm, double *t, double *p, Tatlas *a)
 Samples a random point on the atlas. More...
 
double AtlasVolume (Tparameters *pr, boolean collisionFree, Tatlas *a)
 Approximates the volume of the manifold. More...
 
void SaveChartCenters (Tparameters *pr, char *fname, boolean inside, boolean saveWithDummies, boolean middlePoint, Tatlas *a)
 Stores the centers of the charts. More...
 
void SaveSingularCharts (Tparameters *pr, char *fname, boolean saveWithDummies, Tatlas *a)
 Stores the centers of the singular charts. More...
 
void PlotAtlas (char *fname, int argc, char **arg, Tparameters *pr, FILE *fcost, unsigned int xID, unsigned int yID, unsigned int zID, Tatlas *a)
 Pots a projection of an atlas. More...
 
void TriangulateAtlas (char *fname, int argc, char **arg, Tparameters *pr, FILE *fcost, unsigned int xID, unsigned int yID, unsigned int zID, Tatlas *a)
 Pots the triangular mesh defined by the an atlas. More...
 
unsigned int AtlasMemSize (Tatlas *a)
 Memory used by a given atlas. More...
 
void SaveAtlas (Tparameters *pr, Tfilename *fname, Tatlas *a)
 Stores the atlas information on a file. More...
 
void LoadAtlas (Tparameters *pr, Tfilename *fname, TAtlasBase *w, Tatlas *a)
 Defines an atlas from the information on a file. More...
 
void DeleteAtlas (Tatlas *a)
 Destructor. More...
 

Macro Definition Documentation

#define ATLAS_VERBOSE   0

Vebosity of the atlas operations. If set to 0 only minimalistic information is printed.

Definition at line 29 of file atlas.h.

Referenced by AtlasAStar(), AtlasGBF(), main(), and ReconstructAtlasPath().

#define PLOT_ATLAS_IN_COLORS   1

When set to 1 charts are plotted in different colors according to whether they are singular (green), open (red), singular open (yellow), close (blue), etc. Otherwise all charts are plotted in blue.

This is only used if the user does not provide a cost map associated with the atlas.

Definition at line 39 of file atlas.h.

Referenced by PlotAtlas().

#define PROJECT_WITH_PLANE   1

Select between the two types of projection when interpolating between configurations. In the first type of projection a Newton (or Jacobian pseudo-inverse) method is used. This method does not warrantee any separation between the samples when projected. When we force the projection to lie in a plane orthogonal to the interpolation direction, projected samples are at least as separated as the interpolation steps.

Right now this is only used when approximatting the geodesic distance between samples that is part of the A* search (used to determine if there is a collision between the centers of two charts: collision = INF distance).

Definition at line 53 of file atlas.h.

#define INIT_NUM_CHARTS   100

Initial number of charts of an atlas. It will be expanded as needed.

Definition at line 60 of file atlas.h.

Referenced by AddChart2AtlasRRT(), InitAtlas(), and InitAtlasRRT().

#define INIT_NUM_BIFURCATIONS   100

Initial number of bifurcations of an atlas. It will be expanded as needed.

Definition at line 66 of file atlas.h.

Referenced by InitAtlas().

#define USE_ATLAS_TREE   0

Set to 1 to use a binary tree to search for nighbouring charts. In general should be 1. For very small problems 0 could be better but the difference in execution time will be minor.

Definition at line 78 of file atlas.h.

Referenced by main().

#define SAMPLING_RADIUS_REDUCTION_FACTOR   0.9

Ratio at which the sampling radius is reduced when sampling fails.

Definition at line 85 of file atlas.h.

Referenced by ExtendAtlasFromPoint(), and NewChartFromPoint().

#define SIMPLE_BORDER   1

With this flag on only a rough approximation of the domain borders is performed. It is more efficient.

Definition at line 92 of file atlas.h.

#define INIT_NUM_CHARTS_IN_ATLAS_HEAP   100

Initial number of elements in the heap of local chart used int AtlasGBF. This is extended as needed.

Definition at line 102 of file atlas.h.

Referenced by AtlasAStar(), AtlasGBF(), and BuildAtlasFromPoint().

Function Documentation

void InitAtlasHeapElement ( unsigned int  mID,
double  c,
double  beta,
TAtlasHeapElement he 
)

Defines a new atlas heap element.

Parameters
mIDChart identifier.
cCost (includes the cost from origin plus the heuristic to goal, if any).
betaThe penalization factor (>1).
heThe element to initialize.

Definition at line 2600 of file atlas.c.

References TAtlasHeapElement::beta, TAtlasHeapElement::chartID, TAtlasHeapElement::cost, Error(), and TAtlasHeapElement::nPenalized.

Referenced by AtlasAStar(), AtlasGBF(), and BuildAtlasFromPoint().

void CopyAtlasHeapElement ( void *  he1,
void *  he2 
)

Defines a new atlas heap element from another. We use void pointers to use the generic heap implementation.

Parameters
he1The heap element to define.
he2The heap element from where to copy.

Definition at line 2610 of file atlas.c.

Referenced by AtlasAStar(), AtlasGBF(), and BuildAtlasFromPoint().

unsigned int GetAtlasHeapElementID ( TAtlasHeapElement he)

Gets the chart identifer stored in an atlas heap element.

Parameters
heThe heap element to query.
Returns
The chart identifier stored in the elemement.

Definition at line 2618 of file atlas.c.

References TAtlasHeapElement::chartID.

Referenced by AtlasAStar(), AtlasGBF(), and BuildAtlasFromPoint().

double GetAtlasHeapElementCost ( TAtlasHeapElement he)

Gets the cost stored in an atlas heap element.

Parameters
heThe heap element to query.
Returns
The cost stored in the elemement.

Definition at line 2623 of file atlas.c.

References TAtlasHeapElement::cost.

boolean LessThanAtlasHeapElement ( void *  he1,
void *  he2,
void *  userData 
)

Identifies atlas heap element with lower cost.

Parameters
he1The first heap element to compare.
he2The second heap element to compare.
userDataNot used. Included for compatibility with the generic heap implementation.
Returns
TRUE if he1 is cheaper than he2.
See Also
Theap.

Definition at line 2628 of file atlas.c.

Referenced by AtlasAStar(), AtlasGBF(), and BuildAtlasFromPoint().

void PenalizeAtlasHeapElement ( TAtlasHeapElement he)

Penalizes the cost of the atlas heap element. Every time this function is called we increase the counter in the heap element and, thus, we penalize the cost accordingly. The idea is that charts where it is difficult to sample should be penalized to deviate the search to neighbouring charts with larger external border.

Parameters
heThe atlas heap element to updated.

Definition at line 2643 of file atlas.c.

References TAtlasHeapElement::nPenalized.

Referenced by AtlasGBF().

void DeleteAtlasHeapElement ( void *  he)

Releases the memory alloated in the atlas heap element, if any.

Parameters
heThe heap element to remove.

Definition at line 2648 of file atlas.c.

Referenced by AtlasAStar(), AtlasGBF(), and BuildAtlasFromPoint().

void InitAtlas ( Tparameters pr,
boolean  parallel,
boolean  simpleChart,
unsigned int  k,
double  e,
double  ce,
double  r,
TAtlasBase w,
Tatlas a 
)

Initializes an empty atlas.

Parameters
prThe set of parameters.
parallelTRUE if the atlas might be generated/processed in parallel. This parameter is used if OpenMP is available.
simpleChartTRUE if the atlas is to be defined using simple charts.
kThe dimension of the manifold
eMaximum error between charts and the manifold.
ceMaximum curvature error
rThe radius of validity of the local parametrization.
wThe world on which the atlas is to be defined.
aThe atlas to initialize.

Definition at line 2676 of file atlas.c.

References Tatlas::ambient, Tatlas::bifurcation, Tatlas::ce, Tatlas::charts, CS_WD_GENERATE_SIMP_INITIAL_BOX, CS_WD_GET_SIMP_JACOBIAN, CS_WD_INIT_CD, Tatlas::currentChart, Tatlas::e, FALSE, GetJacobianSize(), Tatlas::H, INIT_NUM_BIFURCATIONS, INIT_NUM_CHARTS, Tatlas::J, Tatlas::k, Tatlas::m, Tatlas::maxCharts, Tatlas::mBifurcations, Tatlas::n, Tatlas::nBifurcations, Tatlas::ncJ, Tatlas::nCores, NEW, Tatlas::npBifurcations, Tatlas::npCharts, Tatlas::nrJ, Tatlas::parallel, Tatlas::r, SetAtlasTopology(), Tatlas::simpleChart, and Tatlas::w.

Referenced by InitAtlasFromPoint().

unsigned int AddChart2Atlas ( Tparameters pr,
double *  ps,
unsigned int  parentID,
boolean singularity,
Tatlas a 
)

Defines a chart on the given point and adds it to the atlas updating the neighbouring relations and checking for the presence of singularities in between the centers of the charts dectected as neighbours, if any.

Note that we could actually build an atlas by sampling random points on the manifold and using this function to create new charts until the atlas is compleated. However this estrategy is not guaranteed be be succesful and could produce atlas with charts for very different sizes. This is why higher dimensional continuation typically uses more principled ways to generate new charts from the existing ones (see BuildAtlasFromPoint).

Parameters
prThe set of parameters.
psThe point (on the manifold) where to define the new chart. CAUTION. This point is assumed to be in the simplified system. (no in the original system as the point given in InitAtlasFromPoint). This kind of points typically result from any manipulation of an existing atlas (RandomPointInAtlas, Chart2Manifold, etc).
parentIDIdentifier of the parent chart. NO_UINT if unkown.
singularity[Output] TRUE if the new chart of one of its neighbours is (almost) on a singularity at the given numerical accuracy (see the epsilon parameter). In this case bifurcations might be undetected.
aThe atlas where to add the new chart.
Returns
Identifier of the new chart or NO_UINT if the chart could not be created.

Definition at line 2836 of file atlas.c.

References Tatlas::ambient, Tatlas::ce, Tatlas::charts, Tatlas::currentChart, Tatlas::e, FALSE, InitChart(), Tatlas::J, Tatlas::k, Tatlas::m, Tatlas::maxCharts, MEM_DUP, NEW, NO_UINT, PostProcessNewCharts(), Tatlas::r, Tatlas::simpleChart, Tatlas::tp, TRUE, and Tatlas::w.

Referenced by MinimizeOnAtlas().

unsigned int AddTrustedChart2Atlas ( Tparameters pr,
double *  ps,
unsigned int  parentID,
boolean singularity,
Tatlas a 
)

Version of AddChart2Atlas for atlas in AtlasRRT (construction of a RRT from an atlas and the reverse).

In this special version, we add the chart to the atlas and we determine its neighbours but without any concern about detecting singularities nor enforcing intersection with the parent chart. Moreover, we trust that the input point is in the manifold and that its is collision free (these conditions are checked when generating the point).

Parameters
prThe set of parameters.
psThe point (on the manifold) where to define the new chart. CAUTION. This point is assumed to be in the simplified system. (no in the original system as the point given in InitAtlasFromPoint). This kind of points typically result from any manipulation of an existing atlas (RandomPointInAtlas, Chart2Manifold, etc).
parentIDIdentifier of the parent chart. NO_UINT if unkown.
singularity[Output] TRUE if the new chart of one of its neighbours is (almost) on a singularity at the given numerical accuracy (see the epsilon parameter). In this case bifurcations might be undetected.
aThe atlas where to add the new chart.
Returns
Identifier of the new chart or NO_UINT if the chart could not be created.

Definition at line 2870 of file atlas.c.

References Tatlas::ambient, Tatlas::ce, Tatlas::charts, CT_DETECT_BIFURCATIONS, Tatlas::currentChart, DetectBifurcation(), Tatlas::e, FALSE, GetParameter(), InitTrustedChart(), Tatlas::J, Tatlas::k, Tatlas::m, Tatlas::maxCharts, MEM_DUP, Tatlas::n, NEW, NO_UINT, PostProcessNewCharts(), Tatlas::r, Tatlas::simpleChart, Tatlas::tp, and Tatlas::w.

Referenced by AddChart2AtlasRRT().

boolean InitAtlasFromPoint ( Tparameters pr,
boolean  parallel,
boolean  simpleChart,
double *  p,
TAtlasBase w,
Tatlas a 
)

Initializes an atlas defining a local chart from a given point. If the given point is not on the manifold nothing is added to the atlas.

Parameters
prThe set of parameters.
parallelTRUE if the atlas might be generated/processed in parallel. Only possible if OpenMP is available.
simpleChartTRUE if the atlas is to be defined using simple charts.
pA point on the manifold. A sample including only system variables in the original system of equations (not the simplified one).
wThe base type with the equations (and possibly obstacles).
aThe atlas to create.
Returns
TRUE if the atlas can be actually created. Errors are triggered if the manifold is singular in the given point.

Definition at line 2742 of file atlas.c.

References Tatlas::ambient, Tatlas::ce, ChangeParameter(), Tatlas::charts, CS_WD_GENERATE_SIMPLIFIED_POINT, CS_WD_ORIGINAL_IN_COLLISION, CS_WD_REGENERATE_SOLUTION_POINT, CT_CE, CT_E, CT_EPSILON, CT_N_DOF, CT_R, Tatlas::currentChart, Tatlas::e, Error(), EvaluateTransposedJacobianInVector(), FindRank(), FrontierChart(), GetParameter(), InitAtlas(), InitBTree(), InitChart(), Tatlas::J, Tatlas::k, Tatlas::m, Tatlas::ncJ, NEW, Tatlas::npCharts, Tatlas::nrJ, Tatlas::r, Tatlas::simpleChart, Tatlas::tp, Tatlas::w, and Warning().

Referenced by BuildAtlasFromPoint(), InitAtlasRRT(), main(), and MinimizeOnAtlas().

void BuildAtlasFromPoint ( Tparameters pr,
double *  p,
boolean  simpleChart,
TAtlasBase w,
Tatlas a 
)

Defines an atlas from a given point and covering the whole connected component including this point.

This procedure can be quite time demanding for high-dimensional manifolds. For these cases, we implement a parallel version (via MP). However, currently the parallel version can not deal with obstacles (i.e., it can only run on cuik files but not on world files).

Todo:
Make collision detection re-entrant so that we can use openMP with it This basically means to have as many collision detection structures as prallel threads.
Parameters
prThe set of parameters.
pA point on the manifold. A sample including only system variables. This is a point in the original system, with values for all the variables.
simpleChartTRUE if the atlas is to be defined using simple charts.
wThe base type with the equations (and possibly obstacles).
aThe atlas to create.

Definition at line 2919 of file atlas.c.

References AddElement2Heap(), BoundaryPointFromExternalCorner(), Tatlas::charts, CopyAtlasHeapElement(), CT_ATLASGBF_BETA, CT_MAX_CHARTS, CT_N_DOF, Tatlas::currentChart, DeleteAtlasHeapElement(), DeleteHeap(), DistanceTopology(), Error(), ExtendAtlasFromPoint(), ExtractMinElement(), FALSE, GetAtlasHeapElementID(), GetChartCenter(), GetParameter(), HeapEmpty(), INIT_NUM_CHARTS_IN_ATLAS_HEAP, InitAtlasFromPoint(), InitAtlasHeapElement(), InitAtlasStatistics(), InitHeap(), Tatlas::k, LessThanAtlasHeapElement(), Tatlas::m, Tatlas::maxCharts, MEM_DUP, Tatlas::nCores, NEW, NewBoundaryAttempt(), NewChartFromPoint(), NewNotInBoundary(), NO_UINT, OpenChart(), Tatlas::parallel, PostProcessNewCharts(), PrintAtlasStatistics(), randomDouble(), RANDOMNESS, Tatlas::tp, TRUE, and WrongCorner().

Referenced by main().

boolean MinimizeOnAtlas ( Tparameters pr,
char *  fname,
double *  p,
TAtlasBase w,
double(*)(Tparameters *, boolean, double *, void *)  costF,
void *  costData,
Tatlas a 
)

Local minimization of a cost function defined on the atlas. The minimization is done using gradient descent. The gradient is numerically evaluated. The presence of bifurcation si taken into account (if the appropiate parameteres are set, see CT_DETECT_BIFURCATIONS) and, thus, different minima can be determined. The paths to the minima are stored in separate files.

Parameters
prThe set of parameters.
fnameBase name used for the output paths.
pThe initial point.
wThe base type with the equations (and possibly obstacles).
costFThe cost function.
costDataThe last parameter of the cost funtion.
aThe atlas to be generated during minimization
Returns
TRUE if the minimization detected no issue. FALSE if the minimization stepped on a singularity and it could not be determined if it corresponds to a bifurcation. In this case the minimization process might be incomplete.

Definition at line 3106 of file atlas.c.

References AddChart2Atlas(), AddSample2Samples(), ChangeParameter(), Chart2Manifold(), Tatlas::charts, CS_WD_GET_SYSTEM_VARS, CS_WD_REGENERATE_ORIGINAL_POINT, CT_DELTA, CT_DETECT_BIFURCATIONS, CT_EPSILON, CT_N_DOF, Tatlas::currentChart, Error(), FALSE, GetChartCenter(), GetParameter(), INF, InitAtlasFromPoint(), InitSamples(), Tatlas::J, Tatlas::m, MEM_DUP, TMinTrace::nc, NEW, Norm(), TMinTrace::ns, SaveSamples(), SaveSamplesN(), ScaleVector(), TMinTrace::st, TRUE, Tatlas::w, and Warning().

Referenced by main().

double GetAtlasRadius ( Tatlas a)

Returns the radius used to initialize the charts.

Parameters
aThe atlas to query.
Returns
The radius used to initialize all radius.

Definition at line 3311 of file atlas.c.

References Tatlas::r.

double GetAtlasError ( Tatlas a)

Returns the maximum error from the chart to the manifold.

Parameters
aThe atlas to query.
Returns
The maximum error used to initialize all the chart in the atlas.

Definition at line 3316 of file atlas.c.

References Tatlas::e.

double GetAtlasErrorCurv ( Tatlas a)

Returns the maximum curvature error for all the charts, i.e., the maximum angle between neighbouring charts.

Parameters
aThe atlas to query.
Returns
The maximum curvature error used to initialize all the chart in the atlas.

Definition at line 3321 of file atlas.c.

References Tatlas::ce.

TAtlasBase* GetAtlasWorld ( Tatlas a)

Returns the workd structure defining the manifold to be described by the atlas.

Parameters
aThe atlas to query.
Returns
The world structure with the manifold description.

Definition at line 3326 of file atlas.c.

References Tatlas::w.

unsigned int GetAtlasAmbientDim ( Tatlas a)

Returns the number of variables of the problem.

Parameters
aThe atlas to query.
Returns
The dimensionality of the ambient space.

Definition at line 3331 of file atlas.c.

References Tatlas::m.

Referenced by main().

unsigned int GetAtlasManifoldDim ( Tatlas a)

Returns the dimensionality of the manifold described by the atlas.

Parameters
aThe atlas to query.
Returns
The dimensionality of the manifold.

Definition at line 3336 of file atlas.c.

References Tatlas::k.

Referenced by main().

boolean Newton2ManifoldPlane ( Tparameters pr,
double *  point,
double *  vector,
double *  pInit,
double *  p,
Tatlas a 
)

Determines a point in the intersection of a given plane and the manifold. This is an auxiliary function of FindPointInOtherBranch but it is also used when defining an (Atlas)RRT (to incrementally extend a branch of the tree toward a given point).

The point is searched using a Newton procedure where the syste of equations defining the manifold is extended with the equation defining the given plane.

Parameters
prThe set of parameters.
pointA passing poing for the plane.
vectorThe vector defining the plane.
pInitInitial estimation for the solution. If NULL, point is used.
pThe solution point, if any. The space for this point must be allocated externally.
aThe atlas.
Returns
TRUE if the process converged to a solution.

Definition at line 2250 of file atlas.c.

References ArrayPi2Pi(), CS_WD_EVALUATE_SIMP_EQUATIONS, CT_EPSILON, CT_MAX_NEWTON_ITERATIONS, DeleteNewton(), Error(), EvaluateJacobianInVector(), FALSE, GetNewtonMatrixBuffer(), GetNewtonRHBuffer(), GetParameter(), InitNewton(), Tatlas::J, Tatlas::m, Tatlas::n, Tatlas::ncJ, NewtonStep(), Norm(), Tatlas::nrJ, SetRow(), Tatlas::tp, TRUE, and Tatlas::w.

Referenced by FindPointInOtherBranch(), and GeodesicDistance().

boolean ExtendAtlasFromPoint ( Tparameters pr,
unsigned int  parentID,
double *  t,
TAtlasStatistics st,
Tatlas a 
)

Tries to expand the atlas from a given point defined on the chart tangent space. This point is assumed to be on the border of the atlas (on the ball defined in the given chart and not in any other ball from other charts).

If the new chart can not be properly defined, the we try to try to generate it closer to the center of the parent chart. When very close to the parent chart center we are supposed to be able to generate a new chart without any issues.

Parameters
prThe set of parameters.
parentIDThe identifier of the chart where the point is defined. This is the parent of the new chart, if we manage to create it.
tThe point in the tangent space to use to expand the atlas..
stTo collect statitic information on the atlas extension process.
aThe atlas to extend.
Returns
TRUE if the altas could be actually extended.

Definition at line 1015 of file atlas.c.

References AddBorderConstraint(), Tatlas::ambient, Tatlas::ce, Chart2Manifold(), ChartIsFrontier(), Tatlas::charts, CloseCharts(), CollisionChart(), CS_WD_IN_COLLISION, CS_WD_SIMP_INEQUALITIES_HOLD, CT_EPSILON, Tatlas::currentChart, DeleteChart(), Tatlas::e, Error(), FALSE, GetParameter(), InitChart(), Tatlas::J, Tatlas::k, Tatlas::m, Tatlas::maxCharts, MEM_DUP, MIN_SAMPLING_RADIUS, NEW, NewAtlasExtension(), NewDecompositionError(), NewFarFromParent(), NewGoodExtension(), NewImpossible(), NewLargeError(), NewNonRegularPoint(), NewNotInManifold(), NewRadiousChange(), NewSingularImpossible(), NO_UINT, Norm(), PointInBoxTopology(), PostProcessNewCharts(), Tatlas::r, SAMPLING_RADIUS_REDUCTION_FACTOR, ScaleVector(), Tatlas::simpleChart, SingularChart(), Tatlas::tp, TRUE, and Tatlas::w.

Referenced by AtlasAStar(), and BuildAtlasFromPoint().

boolean ExtendAtlasTowardPoint ( Tparameters pr,
unsigned int  parentID,
double *  t,
boolean  collisionStops,
TAtlasStatistics st,
Tatlas a 
)

This is like ExtendAtlasFromPoint but the process is started close to the center of the parent chart and the new chart is progressively generated towards the targed point as we have success in generating the new chart. In this process, we could even go further that the targed point.

Another difference with ExtendAtlasFromPoint is that here we take into account obstacles.

Parameters
prThe set of parameters.
parentIDThe identifier of the chart where the point is defined. This is the parent of the new chart, if we manage to create it.
tThe point in the tangent space to use to expand the atlas.
collisionStopsSet to TRUE if the extension has to be stopped when we detect a collision with an obstacle. Otherwise, the extension proceed while the error is moderate.
stTo collect statitic information on the atlas extension process.
aThe atlas to extend.
Returns
TRUE if we could actually extend the atlas.

Definition at line 1246 of file atlas.c.

References Tatlas::ambient, Tatlas::ce, Chart2Manifold(), Tatlas::charts, CloseCharts(), CollisionChart(), CT_DELTA, CT_EPSILON, Tatlas::currentChart, DeleteChart(), Tatlas::e, Error(), FALSE, GetChartCenter(), GetChartRadius(), GetParameter(), InitChart(), Tatlas::J, Tatlas::k, Tatlas::m, Tatlas::maxCharts, MEM_DUP, NEW, NewAtlasExtension(), NewDecompositionError(), NewFarFromParent(), NewGoodExtension(), NewImpossible(), NewInCollision(), NewLargeError(), NewNonRegularPoint(), NewNotInManifold(), NewRadiousChange(), NewSingularImpossible(), Norm(), PostProcessNewCharts(), Tatlas::r, ScaleVector(), Tatlas::simpleChart, SingularChart(), Tatlas::tp, TRUE, and Tatlas::w.

Referenced by AtlasGBF().

boolean AtlasAStar ( Tparameters pr,
double *  p,
double *  time,
double *  pl,
unsigned int *  ns,
double ***  path,
Tatlas a 
)

Expands the atlas with the purpose to reach a given configuration using a A* search strategy.

Parameters
prThe set of parameters.
pThe point to reach.
timeThe actual time used in planning.
plPath lenght.
nsThe number of steps of the solution path. Zero if no solution exists.
pathSequence of ns points to go from the target point to the center of the first chart in the atlas.
aThe atlas to extend.
Returns
TRUE if p could be actually reached.

< Predecesor. Previou node in the best path to the start node.

Definition at line 3341 of file atlas.c.

References AddElement2Heap(), ATLAS_VERBOSE, BoundaryPointFromExternalCorner(), ChartNeighbourID(), ChartNumNeighbours(), Tatlas::charts, CopyAtlasHeapElement(), TAStarInfo::cost, CS_WD_ERROR_IN_SIMP_EQUATIONS, CS_WD_GENERATE_SIMPLIFIED_POINT, CS_WD_GET_SYSTEM_VARS, CS_WD_ORIGINAL_IN_COLLISION, CS_WD_REGENERATE_SOLUTION_POINT, CT_ATLASGBF_BETA, CT_EPSILON, CT_MAX_CHARTS, CT_MAX_PLANNING_TIME, Tatlas::currentChart, DeleteAtlasHeapElement(), DeleteHeap(), DeleteStatistics(), DistanceOnChart(), DistanceTopology(), Error(), ExpandibleChart(), ExtendAtlasFromPoint(), ExtractMinElement(), FALSE, GeodesicDistance(), GetAtlasHeapElementID(), GetChartCenter(), GetElapsedTime(), GetHeapElement(), GetParameter(), HeapEmpty(), HeapSize(), TAStarInfo::heuristic, INF, INIT_NUM_CHARTS_IN_ATLAS_HEAP, InitAtlasHeapElement(), InitAtlasStatistics(), InitHeap(), InitStatistics(), Tatlas::J, Tatlas::k, LessThanAtlasHeapElement(), Tatlas::m, Tatlas::maxCharts, MEM_DUP, MEM_EXPAND, Tatlas::n, Tatlas::nCores, NEW, NewBoundaryAttempt(), NewChartFromPoint(), NewNotInBoundary(), NO_UINT, ON_CUIKSYSTEM, Tatlas::parallel, PointOnChart(), PostProcessNewCharts(), PrintAtlasStatistics(), ReconstructAtlasPath(), TAStarInfo::status, Tatlas::tp, TRUE, Tatlas::w, Warning(), and WrongCorner().

Referenced by main().

boolean AtlasGBF ( Tparameters pr,
double *  p,
double *  time,
double *  pl,
unsigned int *  ns,
double ***  path,
Tatlas a 
)

Expands the atlas with the purpose to reach a given configuration using a Greedy Best First search strategy.

Parameters
prThe set of parameters.
pThe point to reach.
timeThe actual time used in planning.
plPath lenght.
nsThe number of steps of the solution path. Zero if no solution exists.
pathSequence of ns points to go from the target point to the center of the first chart in the atlas.
aThe atlas to extend.
Returns
TRUE if p could be actually reached.

Definition at line 3677 of file atlas.c.

References AddElement2Heap(), ATLAS_VERBOSE, Tatlas::charts, CopyAtlasHeapElement(), CS_WD_ERROR_IN_SIMP_EQUATIONS, CS_WD_GENERATE_SIMPLIFIED_POINT, CS_WD_GET_SYSTEM_VARS, CS_WD_ORIGINAL_IN_COLLISION, CS_WD_REGENERATE_SOLUTION_POINT, CT_ATLASGBF_BETA, CT_EPSILON, CT_MAX_CHARTS, CT_MAX_PLANNING_TIME, Tatlas::currentChart, DeleteAtlasHeapElement(), DeleteHeap(), DeleteStatistics(), DistanceOnChart(), DistanceTopology(), Error(), ExpandibleChart(), ExtendAtlasTowardPoint(), ExtractMinElement(), FALSE, GetAtlasHeapElementID(), GetChartCenter(), GetElapsedTime(), GetParameter(), HeapEmpty(), INF, INIT_NUM_CHARTS_IN_ATLAS_HEAP, InitAtlasHeapElement(), InitAtlasStatistics(), InitHeap(), InitStatistics(), Tatlas::J, Tatlas::k, LessThanAtlasHeapElement(), Tatlas::m, Tatlas::maxCharts, MEM_EXPAND, Tatlas::n, Tatlas::nCores, NEW, NewBoundaryAttempt(), NewNotInBoundary(), NO_UINT, PenalizeAtlasHeapElement(), PointOnChart(), PrintAtlasStatistics(), RandomPointOnBoundary(), ReconstructAtlasPath(), Tatlas::tp, TRUE, Tatlas::w, and Warning().

Referenced by main().

void AtlasTree ( Tparameters pr,
unsigned int  nNodes,
Tatlas a 
)

Tries to define a tree covereing the manifold. \

Parameters
prThe set of parameters.
nNodesThe number of charts (or nodes) of the tree.
aThe atlas to defined. Must have a chart from where to start.
unsigned int GetAtlasNumCharts ( Tatlas a)

Returns the number of charts in the atlas.

Parameters
aThe atlas to query.
Returns
The number of charts in the atlas.

Definition at line 3931 of file atlas.c.

References Tatlas::currentChart.

Referenced by AddChart2AtlasRRT(), InitAtlasRRT(), main(), and RandomPointInAtlasTree().

Tchart* GetAtlasChart ( unsigned int  id,
Tatlas a 
)

Returns a pointer to one of the charts of the atlas. Care should be taken not to manipulate the chart (use it for query only).

Parameters
idThe identifier of the chart.
aThe atlas to query.
Returns
A pointer to the requested chart of NULL if no charts with the given identifier exists.

Definition at line 3936 of file atlas.c.

References Tatlas::charts.

Referenced by AddBranchToAtlasRRT(), AddChart2AtlasRRT(), AddSample2AtlasRRT(), AtlasRRT(), AtlasRRTSample(), GetRRTNNInNeighbourChart(), NewTemptativeSample(), PlotAtlasRRT(), PlotQrand(), PointTowardRandSample(), PopulateWithSamples(), PrintAtlasRRTStatistics(), and RandomPointInAtlasTree().

boolean RandomPointInAtlas ( Tparameters pr,
double  scale,
double *  w,
unsigned int *  nm,
double *  t,
double *  p,
Tatlas a 
)

Samples a random point on the part of the manifold covered by the atlas with uniform distribution.

We select a chart with uniform distribution and then a point in the tangent space of the selected chart.

Parameters
prThe set of parameters.
scaleGlobal scale for sampling areas of the charts. This is only used for simple charts (see RandomPointInChart).
wWeight for the charts. If NULL uniform distribution is used.
nmIdentifier of the chart used to generate the random point.
tParameters in the chart for the random point.
pThe output point in the ambient space (the point int he chart given by t possibly projected to the manifold).
aThe atlas to sample.
Returns
TRUE if we actually managed to sample a point in the atlas with the given requirements.

Definition at line 3944 of file atlas.c.

References Tatlas::charts, Tatlas::currentChart, randomMax(), RandomPointInChart(), randomWithDistribution(), and Tatlas::tp.

Referenced by RandomPointInAtlasTree().

double AtlasVolume ( Tparameters pr,
boolean  collisionFree,
Tatlas a 
)

Approximates the volume of the manifold parametrized by an atlas.

Parameters
prThe set of parameters.
collisionFreeTRUE if only volume of the collision free of the manifold has to be approximated.
aThe atlas.
Returns
The volume of the manifold.

Definition at line 3966 of file atlas.c.

References Tatlas::charts, ChartVolume(), Tatlas::currentChart, FrontierChart(), Tatlas::J, and Tatlas::tp.

Referenced by main().

void SaveChartCenters ( Tparameters pr,
char *  fname,
boolean  inside,
boolean  saveWithDummies,
boolean  middlePoint,
Tatlas a 
)

Stores the centers of the charts in the form of boxes.

Parameters
prThe set of parameters.
fnameName to use for the two output file.
insideTRUE if only points inside the domain have to be saved (i.e., exclude FrontierCharts)
saveWithDummiesTRUE if the output has to include de dummies. This is necessary only if you plan to unsimplify the output points.
middlePointIf TRUE, we output includes an intermediate point in between each pair of neighbouring charts.
aThe atlas to query.

Definition at line 3987 of file atlas.c.

References Chart2Manifold(), ChartNeighbourID(), ChartNumNeighbours(), Tatlas::charts, CreateFileName(), CS_WD_GET_SYSTEM_VARS, CS_WD_REGENERATE_ORIGINAL_POINT, Tatlas::currentChart, DeleteFileName(), Error(), FrontierChart(), GetChartCenter(), GetFileFullName(), INF, Tatlas::J, Tatlas::k, Tatlas::m, Manifold2Chart(), NEW, NO_UINT, ScaleVector(), SOL_EXT, SOL_WITH_DUMMIES_EXT, Tatlas::tp, and Tatlas::w.

Referenced by main().

void SaveSingularCharts ( Tparameters pr,
char *  fname,
boolean  saveWithDummies,
Tatlas a 
)

Stores the centers of the singular charts in the form of boxes.

Parameters
prThe set of parameters.
fnameName to use for the two output file.
saveWithDummiesTRUE if the output has to include de dummies. This is necessary only if you plan to unsimplify the output points.
aThe atlas to query.

Definition at line 4095 of file atlas.c.

References Tatlas::charts, CreateFileName(), CS_WD_GET_SYSTEM_VARS, CS_WD_REGENERATE_ORIGINAL_POINT, Tatlas::currentChart, DeleteFileName(), Error(), GetChartCenter(), GetFileFullName(), SingularChart(), SOL_EXT, SOL_WITH_DUMMIES_EXT, and Tatlas::w.

Referenced by main().

void PlotAtlas ( char *  fname,
int  argc,
char **  arg,
Tparameters pr,
FILE *  fcost,
unsigned int  xID,
unsigned int  yID,
unsigned int  zID,
Tatlas a 
)

Plots a 3d projection of an atlas defined on a manifold.

Although the ambien space can have arbitrary dimension we project it on 3 dimensions. Only 2D manifolds plots can be properly visualized.

The output plot can be visualized using geomview.

Parameters
fnameFile where to store the plot.
argcNumber arguments given to the program calling this function. This is used to write commnets in the output gcl file so that the plot can be easiy reproduced.
argArguments given to the program calling this function. This is used This is used to write commnets in the output gcl file so that the plot can be easiy reproduced.
prThe set of parameters.
fcostOne cost for each chart in the atlas. Used to associate cost maps with the atlas representing each chart with a different color. Set to NULL if the atlas is to have a uniform color.
xIDThe first ambient dimension where to project (in the original system including system vars and dummies).
yIDThe second ambient dimension where to project (in the original system including system vars and dummies).
zIDThe thrid ambient dimension where to project (in the original system including system vars and dummies).
aThe atlas to plot.
See Also
PlotChart.

Definition at line 4151 of file atlas.c.

References ChartNeighbourID(), ChartNumNeighbours(), Tatlas::charts, Close3dObject(), Close3dObjectNoColor(), ClosePlot3d(), CS_WD_REGENERATE_ORIGINAL_POINT, Tatlas::currentChart, Delete3dObject(), DeleteColor(), ExpandibleChart(), FALSE, FrontierChart(), GetChartCenter(), InitPlot3d(), Tatlas::J, Tatlas::k, NewColor(), NO_UINT, PLOT_AS_POLYGONS, PLOT_ATLAS_IN_COLORS, PlotChart(), PlotVect3d(), Tatlas::simpleChart, SingularChart(), StartNew3dObject(), and Tatlas::w.

Referenced by main(), and PlotAtlasRRT().

void TriangulateAtlas ( char *  fname,
int  argc,
char **  arg,
Tparameters pr,
FILE *  fcost,
unsigned int  xID,
unsigned int  yID,
unsigned int  zID,
Tatlas a 
)

Plots a 3d projection of triangular mesh formed by the atlas centers. This is only well defined for surfaces, i.e., for manifolds of dimension 2.

The difference with PlotAtlas is that here only the chart centers are used and, thus, the charts are not displayed. The use of a triangular mesh on the chart centers results in a smoother visualization. The plot generated with PlotAtlas creates a separate polytope for each atlas and, thus there are "discontinuities" in between charts. The triangular mesh is fully continuous and can be properly smoothed and shaded.

Although the ambien space can have arbitrary dimension we project it on 3 dimensions.

The output plot can be visualized using geomview.

Parameters
fnameFile where to store the plot.
argcNumber arguments given to the program calling this function. This is used to write commnets in the output gcl file so that the plot can be easiy reproduced.
argArguments given to the program calling this function. This is used This is used to write commnets in the output gcl file so that the plot can be easiy reproduced.
prThe set of parameters.
fcostOne cost for each chart in the atlas. Used to associate cost maps with the atlas representing each chart with a different color. Set to NULL if the atlas is to have a uniform color.
xIDThe first ambient dimension where to project (in the original system including system vars and dummies).
yIDThe second ambient dimension where to project (in the original system including system vars and dummies).
zIDThe thrid ambient dimension where to project (in the original system including system vars and dummies).
aThe atlas to plot.
See Also
PlotAtlas.

Definition at line 4314 of file atlas.c.

References Tatlas::charts, Close3dObject(), ClosePlot3d(), CostColor(), CrossProduct(), CrossTopologyBorder(), CS_WD_GET_VAR_TOPOLOGY, CS_WD_REGENERATE_ORIGINAL_POINT, CT_CUT_X, CT_CUT_Y, CT_CUT_Z, CT_REPRESENTATION, Tatlas::currentChart, DifferenceVector(), DotProduct(), Error(), FALSE, GetChartCenter(), GetChartNeighboursFromVertices(), GetParameter(), InitPlot3d(), Tatlas::k, M_2PI, MEM_DUP, NeighbouringTriangles(), NEW, NewColor(), NO_UINT, Plot3dObject(), Plot3dObjectWithColors(), REP_JOINTS, SameTriangle(), ScaleVector(), StartNew3dObject(), TOPOLOGY_R, TRUE, and Tatlas::w.

Referenced by main().

unsigned int AtlasMemSize ( Tatlas a)

Returns the approximated memory used (in bytes) by a given atlas.

Parameters
aThe atlas.
Returns
The size of the atlas in bytes.

Definition at line 3831 of file atlas.c.

References ChartMemSize(), Tatlas::charts, and Tatlas::currentChart.

Referenced by AtlasRRTMemSize(), and main().

void SaveAtlas ( Tparameters pr,
Tfilename fname,
Tatlas a 
)

Stores all the information in the atlas in a file.

Todo:
Save atlas (and charts) in binary format.
Parameters
prThe set of parameters.
fnameName of the file where to store the atlas.
aThe atlas to store.

Definition at line 3842 of file atlas.c.

References Tatlas::ce, Tatlas::charts, Tatlas::currentChart, Tatlas::e, Error(), GetFileFullName(), Tatlas::k, Tatlas::m, Tatlas::maxCharts, Tatlas::n, Tatlas::r, SaveBifurcations(), SaveChart(), and Tatlas::simpleChart.

Referenced by main(), and SaveAtlasRRT().

void LoadAtlas ( Tparameters pr,
Tfilename fname,
TAtlasBase w,
Tatlas a 
)

Construct an atlas from the information previously stored in a file by SaveAtlas.

Parameters
prThe set of parameters.
fnameThe name of the file from where to get the information.
wThe world with the equations defining the manifold where the atlas is defined. This is not stored in the file and must be provided externally.
aThe atlas to define.

Definition at line 3869 of file atlas.c.

References AddChart2Btree(), Tatlas::ambient, Tatlas::ce, Tatlas::charts, CS_WD_GENERATE_SIMP_INITIAL_BOX, CS_WD_GET_SIMP_JACOBIAN, Tatlas::currentChart, Tatlas::e, Error(), FALSE, GetFileFullName(), GetJacobianSize(), Tatlas::H, InitBTree(), Tatlas::J, Tatlas::k, LoadBifurcations(), LoadChart(), Tatlas::m, Tatlas::maxCharts, Tatlas::n, Tatlas::ncJ, Tatlas::nCores, NEW, Tatlas::nrJ, Tatlas::parallel, Tatlas::r, SetAtlasTopology(), Tatlas::simpleChart, Tatlas::tp, and Tatlas::w.

Referenced by LoadAtlasRRT(), and main().

void DeleteAtlas ( Tatlas a)

Deletes the information stored in the atlas.

Parameters
aThe atlas to delete.

Definition at line 4561 of file atlas.c.

References Tatlas::ambient, Tatlas::charts, Tatlas::currentChart, DeleteBifurcations(), DeleteBox(), DeleteBTree(), DeleteChart(), DeleteHessian(), DeleteJacobian(), Tatlas::H, Tatlas::J, and Tatlas::tp.

Referenced by DeleteAtlasRRT(), and main().