(*******************************************************************
This file was generated automatically by the Mathematica front end.
It contains Initialization cells from a Notebook file, which
typically will have the same name as this file except ending in
".nb" instead of ".m".

This file is intended to be loaded into the Mathematica kernel using
the package loading commands Get or Needs.  Doing so is equivalent
to using the Evaluate Initialization Cells menu command in the front
end.

DO NOT EDIT THIS FILE.  This entire file is regenerated
automatically each time the parent Notebook file is saved in the
Mathematica front end.  Any changes you make to this file will be
overwritten.
***********************************************************************)









BeginPackage["Simplex`SimplexBasics`"]



Simplex`SimplexBasics::usage="The Simplex`SimplexBasics` package contains \
basic routines for converting points and vectors back and forth between the \
2-d simplex coordinates (x,y,z) with x+y+z=1 and a 2d Cartesian coordinates \
wth the simplex corners at (0,0),(Sqrt[3]/2,1/2), and (0,1)."







Simplex2dCoord::usage="Simplex2dCoord[{x1_,x2_,x3_}] transforms
simplex vector {x1,x2,x3} into the cartesian coordinate vector {y1,y2}
for the simplex defined by corners (0,0), (1/2,sqrt[3]/2), (0,1)."






Simplex2dInverseCoord::usage="Simplex2dInverseCoord[{y1,y2}]
transforms the cartesian coordinate vector {y1,y2} into the simplex
vector {x1,x2,x3} for the simplex defined by corners (0,0),
(1/2,sqrt[3]/2), (0,1)."






Simplex2dVector::usage="Simplex2dVector[{x1_,x2_,x3_},{vx1_,vx2_,vx3_}]
computes the cartesian vector {y1,y2} from point pt={x1,x2,x3} on the
simplex equivalent to the simplex vector {vx1,vx2,vx3}"






DistanceApart::usage="DistanceApart[{x1,y1},{x2,y2}] gives the
distance between the two points."






WithinSimplexBoundary::usage="WithinSimplexBoundary[{x1_,x2_,x3_}] returns 1 \
if
{x1,x2,x3} is in the unit simplex, 0 if otherwise"





WithinSimplexBoundaryNull::usage="WithinSimplexBoundaryNull[{x1_,x2_,x3_}] \
returns 1 if
{x1,x2,x3} is in the unit simplex, Null if otherwise"







WithinSimplexSolution::usage="WithinSimplexSolution[SolutionList] takes a \
list of solutions (e.g. to some set of equations) and returns the first that \
lies within the positive unit simplex of any dimension. "





SimplexUniformRV::usage=
  "SimplexUniformRV[dimension_] uses a broken stick algorithm to draws a \
random point from a uniform distrubtion over the n-dimensional simplex (i.e., \
a uniform Dirichlet distribution)."





Begin["`Private`"]







Simplex2dCoord[{x1_,x2_,x3_}]:=x1*{0,0}+x2*{1,0}+x3*{1/2,Sqrt[3]/2}




Simplex2dInverseCoord[{y1_,y2_}]:={1,0,0}+y1*{-1,1,0}+
    y2*{-1/Sqrt[3],-1/Sqrt[3],2/Sqrt[3]}




DistanceApart[{x1_,y1_},{x2_,y2_}]:=((x1-x2)^2+(y1-y2)^2)^(1/2)



Simplex2dVector[{x1_,x2_,x3_},{vx1_,vx2_,vx3_}]:=
  vx1*(-Simplex2dCoord[{x1,x2,x3}]+Simplex2dCoord[{1,0,0}])/
        DistanceApart[Simplex2dCoord[{x1,x2,x3}],Simplex2dCoord[{1,0,0}]]+
    vx2*(-Simplex2dCoord[{x1,x2,x3}]+Simplex2dCoord[{0,1,0}])/
        DistanceApart[Simplex2dCoord[{x1,x2,x3}],Simplex2dCoord[{0,1,0}]]+
    vx3*(-Simplex2dCoord[{x1,x2,x3}]+Simplex2dCoord[{0,0,1}])/
        DistanceApart[Simplex2dCoord[{x1,x2,x3}],Simplex2dCoord[{0,0,1}]]




WithinSimplexBoundary[{x1_,x2_,x3_}]:=
  If[x1\[GreaterEqual]0&&x2\[GreaterEqual]0&&x3\[GreaterEqual]0,1,0]

WithinSimplexBoundaryNull[{x1_,x2_,x3_}]:=
  If[x1\[GreaterEqual]0&&x2\[GreaterEqual]0&&x3\[GreaterEqual]0,1,Null]



WithinSimplexSolution[solutionList_]:=
  Module[{justValues,testWithinSimplex},
    justValues[solutionl_]:=
      Partition[
        Transpose[Partition[Flatten[Flatten[solutionl,3],1,Rule],2]][[2]],
        Length[solutionl[[1]]]];
    testWithinSimplex[values_]:=
      If[Apply[Plus,values]\[Equal]1&&Apply[Plus,Abs[values]]\[Equal]1,1,0];
    solutionList[[Position[
            Map[testWithinSimplex,justValues[solutionList]],1][[1,1]]]]]




SimplexUniformRV[dimension_] :=
  Partition[Flatten[{0,Sort[Table[Random[],{dimension}]],1}],2,1].{-1,1}




End[]





EndPackage[]