QEIOImpl.mesa
DIRECTORY
Rope,
IO,
Convert,
AlgebraClasses,
Points,
Ints,
BigRats,
Reals,
RatIntervals,
Variables,
DistribPolys,
Polynomials,
AlgebraicNumbers,
ExtensionFields,
FormulaOperators,
Formulas,
QETypes,
QEIO,
Graphs;
QEIOImpl: CEDAR PROGRAM
IMPORTS Rope, IO, Convert, Points, Ints, BigRats, Reals, RatIntervals, Variables, Polynomials, AlgebraicNumbers, ExtensionFields, Formulas, FormulaOperators, Graphs
EXPORTS QEIO
~ BEGIN OPEN AC: AlgebraClasses, PTS: Points, BR: BigRats, RI: RatIntervals, VARS: Variables, DP: DistribPolys, POL: Polynomials, AN: AlgebraicNumbers, EF: ExtensionFields, QFF: Formulas, QET: QETypes, QEIO;
Cylindrical Algebraic Decompositions
ReadCad: PUBLIC PROC [in: IO.STREAM] RETURNS [out: QET.Cad] ~ {
token: Rope.ROPE;
dimension: CARDINAL;
inputVariables, projVariables, minPolyVariable, fieldElementVariable: VARS.VariableSeq;
inputPolynomialRing, projPolynomialRing, minPolyRing, formulaAlgebra: AC.Structure;
inputPolynomials, basis: POL.PolynomialSeq ← NIL;
contents, projection: POL.PolynomialSeq ← NIL;
regions: QET.SignedRegionSeq ← NIL;
cells: QET.CellSeq ← NIL;
cellIndexStructure, inputSignatureStructure: AC.Structure;
inducedCad: QET.Cad ← NIL;
token ← in.GetID[];
IF NOT Rope.Equal[token, "BEGINCAD"] THEN ERROR;
token ← in.GetID[];
IF NOT Rope.Equal[token, "DIMENSION"] THEN ERROR;
dimension ← in.GetCard[];
cellIndexStructure ← PTS.MakePointStructure[Ints.Ints, dimension];
token ← in.GetID[];
IF NOT Rope.Equal[token, "INPUTVARIABLES"] THEN ERROR;
inputVariables ← VARS.ReadVariableSeq[in];
inputPolynomialRing ← POL.MakePolynomialStructure[BR.BigRats, inputVariables];
token ← in.GetID[];
IF NOT Rope.Equal[token, "MINPOLYVARIABLE"] THEN ERROR;
minPolyVariable ← VARS.ReadVariableSeq[in];
minPolyRing ← POL.MakePolynomialStructure[BR.BigRats, minPolyVariable];
token ← in.GetID[];
IF NOT Rope.Equal[token, "FIELDELEMENTVARIABLE"] THEN ERROR;
fieldElementVariable ← VARS.ReadVariableSeq[in];
IF NOT Rope.Equal[fieldElementVariable[1], minPolyVariable[1] ] THEN ERROR; -- quick fix; the entire use of fieldElementVariable should be removed
token ← in.GetID[];
IF Rope.Equal[token, "INPUTPOLYNOMIALS"] THEN {
inputPolynomials ← POL.ReadPolySeq[in, inputPolynomialRing];
inputSignatureStructure ← PTS.MakePointStructure[FormulaOperators.Operators, inputPolynomials.lengthPlus1-1];
token ← in.GetID[];
}
ELSE inputSignatureStructure ← NIL;
IF Rope.Equal[token, "BASIS"] THEN {
basis ← POL.ReadPolySeq[in, inputPolynomialRing];
token ← in.GetID[];
};
projVariables ← VARS.RemoveMainVariable[inputVariables];
projPolynomialRing ← POL.MakePolynomialStructure[BR.BigRats, projVariables];
IF Rope.Equal[token, "CONTENTS"] THEN {
contents ← POL.ReadPolySeq[in, projPolynomialRing];
token ← in.GetID[];
};
IF Rope.Equal[token, "PROJECTION"] THEN {
projection ← POL.ReadPolySeq[in, projPolynomialRing];
token ← in.GetID[];
};
formulaAlgebra ← QFF.MakeFormulaAlgebra[inputPolynomialRing];
IF Rope.Equal[token, "CELLS"] THEN [cells, token] ← ReadCellSeq[in, cellIndexStructure, inputPolynomialRing, inputSignatureStructure, minPolyRing, formulaAlgebra];
IF Rope.Equal[token, "CLUSTERS"] THEN [regions, token] ← ReadSignedRegionSeq[in, cellIndexStructure, inputPolynomialRing, minPolyRing, inputSignatureStructure];
IF Rope.Equal[token, "INDUCEDCAD"] THEN {
inducedCad ← ReadCad[in];
token ← in.GetID[];
};
IF NOT Rope.Equal[token, "ENDCAD"] THEN ERROR;
out ← NEW[QET.CadRec ← [
dimension: dimension,
inputVariables: inputVariables,
minPolyVariable: minPolyVariable,
fieldElementVariable: fieldElementVariable,
inputPolynomials: inputPolynomials,
basis: basis,
contents: contents,
basisProjection: projection,
cells: cells,
clusters: regions,
inducedCad: inducedCad
] ];
FOR i:NAT IN [1..cells.lengthPlus1 - 1] DO
cells[i].cad ← out;
ENDLOOP;
out.abstractGraph ← NEW [Graphs.GraphPrivate ← [
class: cadGraphClass,
rep: out
]];
RETURN;
};
CadFromRope: PUBLIC PROC [in: Rope.ROPE] RETURNS [out: QET.Cad] ~ {
stream: IO.STREAM ← IO.RIS[in];
out ← ReadCad[stream];
};
CadToRope: PUBLIC PROC [in: QET.Cad] RETURNS [out: Rope.ROPE] ~ {
projVariables: VARS.VariableSeq;
out ← "BEGINCAD\n";
out ← Rope.Concat[out, "DIMENSION\n"];
out ← Rope.Cat[out, Convert.RopeFromCard[in.dimension], "\n"];
out ← Rope.Concat[out, "INPUTVARIABLES\n"];
out ← Rope.Cat[out, VARS.VariableSeqToRope[in.inputVariables], "\n"];
out ← Rope.Concat[out, "MINPOLYVARIABLE\n"];
out ← Rope.Cat[out, VARS.VariableSeqToRope[in.minPolyVariable], "\n"];
out ← Rope.Concat[out, "FIELDELEMENTVARIABLE\n"];
out ← Rope.Cat[out, VARS.VariableSeqToRope[in.fieldElementVariable], "\n"];
IF in.inputPolynomials#NIL THEN {
out ← Rope.Concat[out, "INPUTPOLYNOMIALS\n"];
out ← Rope.Cat[out, POL.PolySeqToRope[in.inputPolynomials], "\n"];
};
IF in.basis#NIL THEN {
out ← Rope.Concat[out, "BASIS\n"];
out ← Rope.Cat[out, POL.PolySeqToRope[in.basis], "\n"];
};
projVariables ← VARS.RemoveMainVariable[in.inputVariables];
IF in.contents#NIL THEN {
out ← Rope.Concat[out, "CONTENTS\n"];
out ← Rope.Cat[out, POL.PolySeqToRope[in.contents], "\n"];
};
IF in.basisProjection#NIL THEN {
out ← Rope.Concat[out, "PROJECTION\n"];
out ← Rope.Cat[out, POL.PolySeqToRope[in.basisProjection], "\n"];
};
IF in.cells#NIL THEN {
out ← Rope.Concat[out, "CELLS\n"];
out ← Rope.Cat[out, CellSeqToRope[in.cells], "\n"];
};
IF in.clusters#NIL THEN {
out ← Rope.Concat[out, "CLUSTERS\n"];
out ← Rope.Cat[out, SignedRegionSeqToRope[in.clusters], "\n"];
};
IF in.inducedCad#NIL THEN {
out ← Rope.Concat[out, "INDUCEDCAD\n"];
out ← Rope.Cat[out, CadToRope[in.inducedCad], "\n"];
};
out ← Rope.Concat[out, "ENDCAD\n"];
};
WriteCad: PUBLIC PROC [in: QET.Cad, out: IO.STREAM] ~ {
cadRope: Rope.ROPE ← CadToRope[in];
out.PutRope[cadRope];
};
GetInducedCad: PUBLIC PROC [inCad: QET.Cad, dimension: CARDINAL] RETURNS [outCad: QET.Cad] ~ {
IF inCad.dimension < dimension THEN ERROR;
IF dimension = 0 THEN RETURN[NIL];
outCad ← inCad;
WHILE outCad.dimension > dimension DO
outCad ← outCad.inducedCad;
ENDLOOP;
};
cadVertexClass: PUBLIC Graphs.VertexClass ← Graphs.NewVertexClass[[
Expand: CellExpand
]];
cadGraphClass: PUBLIC Graphs.GraphClass ← Graphs.NewGraphClass[[
EnumerateVertices: CadEnumerateVertices,
Destroy: NIL
]];
CadEnumerateVertices: PUBLIC Graphs.EnumerateVerticesProc ~ {
cad: QET.Cad ← NARROW[graph.rep];
cells: QET.CellSeq ← cad.cells;
FOR i:NAT IN [1..cad.cells.lengthPlus1 - 1] DO
consume.proc[cells[i].abstractVertex, consume.data];
ENDLOOP;
};
CellExpand: PUBLIC Graphs.ExpandProc ~ {
cell: QET.Cell ← NARROW[vertex.rep];
adjacencies: QET.CellSeq ← cell.adjacentCells;
FOR i:NAT IN [1..adjacencies.lengthPlus1 - 1] DO
edge: Graphs.Edge ← [
direction: Undirected,
label: NIL,
otherSide: adjacencies[i].abstractVertex
];
consume.proc[edge, consume.data];
ENDLOOP;
};
Clusters and Cluster Sequences
ReadSignedRegion: PUBLIC PROC [in: IO.STREAM, cellIndexStructure, inputPolynomialRing, minPolyRing, signatureStructure: AC.Structure] RETURNS [out: QET.SignedRegion] ~ {
token: Rope.ROPE ← in.GetID[];
polyRingData: POL.PolynomialRingData ← NARROW[inputPolynomialRing.instanceData];
signature: QFF.Formula;
samplePoint: QET.SamplePoint;
constituentCells: QET.CellIndexSeq;
IF NOT Rope.Equal[token, "SIGNATURE"] THEN ERROR;
signature ← ReadSignature[in, signatureStructure];
token ← in.GetID[];
IF NOT Rope.Equal[token, "SAMPLEPOINT"] THEN ERROR;
samplePoint ← ReadSamplePoint[in, inputPolynomialRing, minPolyRing];
token ← in.GetID[];
IF NOT Rope.Equal[token, "CONSTITUENTCELLS"] THEN ERROR;
constituentCells ← ReadCellIndexSeq[in, cellIndexStructure];
out ← NEW[QET.SignedRegionRec ← [
signature: signature,
samplePoint: samplePoint,
constituentCells: constituentCells
] ];
};
SignedRegionFromRope: PUBLIC PROC [in: Rope.ROPE, cellIndexStructure, inputPolynomialRing, minPolyRing, signatureStructure: AC.Structure] RETURNS [out: QET.SignedRegion] = {
out ← ReadSignedRegion[IO.RIS[in], cellIndexStructure, inputPolynomialRing, minPolyRing, signatureStructure];
};
SignedRegionToRope: PUBLIC PROC [in: QET.SignedRegion] RETURNS [out: Rope.ROPE] ~ {
out ← "CLUSTER\n";
out ← Rope.Concat[out, "SIGNATURE\n"];
out ← Rope.Cat[out, SignatureToRope[in.signature], "\n"];
out ← Rope.Concat[out, "SAMPLEPOINT\n"];
out ← Rope.Concat[out, SamplePointToRope[in.samplePoint]];
out ← Rope.Concat[out, "CONSTITUENTCELLS\n"];
out ← Rope.Cat[out, CellIndexSeqToRope[in.constituentCells], "\n"];
};
WriteSignedRegion: PUBLIC PROC [in: QET.SignedRegion, out: IO.STREAM] = {
out.PutRope[ SignedRegionToRope[in] ]
};
ReadSignedRegionSeq: PUBLIC PROC [in: IO.STREAM, cellIndexStructure, inputPolynomialRing, minPolyRing, signatureStructure: AC.Structure] RETURNS [clusters: QET.SignedRegionSeq, termToken: Rope.ROPE] ~ {
length: NAT ← 0;
token: Rope.ROPE ← in.GetID[];
cluster: QET.SignedRegion;
clusterList, clusterListTail: LIST OF QET.SignedRegion ← NIL;
WHILE Rope.Equal[token, "CLUSTER"] DO
cluster ← ReadSignedRegion[in, cellIndexStructure, inputPolynomialRing, minPolyRing, signatureStructure];
length ← length + 1;
IF clusterList=NIL THEN clusterList ← clusterListTail ←LIST[cluster] ELSE
{ clusterListTail.rest ← LIST[cluster]; clusterListTail ← clusterListTail.rest };
token ← in.GetID[];
ENDLOOP;
clusters ← NEW[QET.SignedRegionSeqRec[length]]; -- assumes length >=1
FOR i:NAT IN [1..length] DO
clusters[i] ← clusterList.first;
clusterList ← clusterList.rest;
ENDLOOP;
RETURN[clusters, token];
};
SignedRegionSeqFromRope: PUBLIC PROC [in: Rope.ROPE, cellIndexStructure, inputPolynomialRing, minPolyRing, signatureStructure: AC.Structure] RETURNS [regions: QET.SignedRegionSeq] ~ {
stream: IO.STREAM ← IO.RIS[in];
termToken: Rope.ROPE;
[regions, termToken ] ← ReadSignedRegionSeq[stream, cellIndexStructure, inputPolynomialRing, minPolyRing, signatureStructure];
};
SignedRegionSeqToRope: PUBLIC PROC [in: QET.SignedRegionSeq] RETURNS [out: Rope.ROPE] ~ {
out ← "";
FOR i:NAT IN [1..in.lengthPlus1 - 1] DO
out ← Rope.Concat[out, SignedRegionToRope[in[i] ] ];
ENDLOOP;
};
WriteSignedRegionSeq: PUBLIC PROC [in: QET.SignedRegionSeq, out: IO.STREAM] ~ {
SCSRope: Rope.ROPE ← SignedRegionSeqToRope[in];
out.PutF["%g", IO.rope[SCSRope] ];
};
Cells and Cell Sequences
ReadCell: PUBLIC PROC [in: IO.STREAM, cellIndexStructure, inputPolynomialRing, inputSignatureStructure, minPolyRing, formulaAlgebra: AC.Structure] RETURNS [out: QET.Cell, termToken: Rope.ROPE] ~ {
token: Rope.ROPE;
index: QET.CellIndex;
inputSignature: QET.Signature;
samplePoint: QET.SamplePoint ← NIL;
definingFormula: QFF.Formula ← NIL;
termChar: CHAR;
token ← in.GetID[];
IF NOT Rope.Equal[token, "INDEX"] THEN ERROR;
index ← ReadCellIndex[in, cellIndexStructure];
token ← in.GetID[];
IF Rope.Equal[token, "INPUTSIGNATURE"] THEN {
inputSignature ← ReadSignature[in, inputSignatureStructure];
token ← in.GetID[];
};
IF NOT Rope.Equal[token, "SAMPLEPOINT"] THEN ERROR;
samplePoint ← ReadSamplePoint[in, inputPolynomialRing, minPolyRing];
token ← in.GetID[];
IF Rope.Equal[token, "DEFININGFORMULA"] THEN {
[definingFormula, termChar] ← QFF.ReadFormula[in, formulaAlgebra];
token ← in.GetID[];
};
IF NOT Rope.Equal[token, "ENDCELL"] THEN ERROR;
token ← in.GetID[];
out ← NEW[QET.CellRec ← [
dimension: CellDimension[index],
index: index,
inputSignature: inputSignature,
samplePoint: samplePoint,
definingFormula: definingFormula
adjacentCells: adjacentCells
]];
out.abstractVertex ← NEW [Graphs.VertexPrivate ← [
class: cadVertexClass,
rep: out
]];
RETURN[out, token];
};
CellFromRope: PUBLIC PROC [in: Rope.ROPE, cellIndexStructure, inputPolynomialRing, inputSignatureStructure, minPolyRing, formulaAlgebra: AC.Structure] RETURNS [out: QET.Cell, termToken: Rope.ROPE] = {
[out, termToken] ← ReadCell[IO.RIS[in], cellIndexStructure, inputPolynomialRing, inputSignatureStructure, minPolyRing, formulaAlgebra];
};
CellToRope: PUBLIC PROC [in: QET.Cell] RETURNS [out: Rope.ROPE] ~ {
out ← "BEGINCELL\n";
out ← Rope.Concat[out, "INDEX\n"];
out ← Rope.Cat[out, CellIndexToRope[in.index], "\n"];
IF in.inputSignature#NIL THEN {
out ← Rope.Concat[out, "INPUTSIGNATURE\n"];
out ← Rope.Cat[out, SignatureToRope[in.inputSignature], "\n"];
};
out ← Rope.Concat[out, "SAMPLEPOINT\n"];
out ← Rope.Concat[out, SamplePointToRope[in.samplePoint]];
IF in.definingFormula#NIL THEN {
out ← Rope.Concat[out, "DEFININGFORMULA\n"];
out ← Rope.Cat[out, QFF.ToRope[in.definingFormula], "\n"];
};
out ← Rope.Concat[out, "ENDCELL\n"];
};
WriteCell: PUBLIC PROC [in: QET.Cell, out: IO.STREAM] = {
out.PutRope[ CellToRope[in] ]
};
LookupCell: PUBLIC PROC [cells: QET.CellSeq, index: QET.CellIndex] RETURNS [QET.Cell] ~ {
FOR i:NAT IN [1..cells.lengthPlus1 - 1] DO
IF CellIndexEqual[cells[i].index, index] THEN RETURN[cells[i] ];
ENDLOOP;
RETURN[NIL];
};
ReadCellSeq: PUBLIC PROC [in: IO.STREAM, cellIndexStructure, inputPolynomialRing, inputSignatureStructure, minPolyRing, formulaAlgebra: AC.Structure] RETURNS [cells: QET.CellSeq, termToken: Rope.ROPE] ~ {
length: NAT ← 0;
token: Rope.ROPE ← in.GetID[];
cell: QET.Cell;
cellList, cellListTail: LIST OF QET.Cell ← NIL;
WHILE Rope.Equal[token, "BEGINCELL"] DO
[cell, token] ← ReadCell[in, cellIndexStructure, inputPolynomialRing, inputSignatureStructure, minPolyRing, formulaAlgebra];
length ← length + 1;
IF cellList=NIL THEN cellList ← cellListTail ←LIST[cell] ELSE
{ cellListTail.rest ← LIST[cell]; cellListTail ← cellListTail.rest };
ENDLOOP;
cells ← NEW[QET.CellSeqRec[length]]; -- assume: "CLUSTERS" seen => length >=1
FOR i:NAT IN [1..length] DO
cells[i] ← cellList.first;
cellList ← cellList.rest;
ENDLOOP;
RETURN[cells, token];
};
CellSeqFromRope: PUBLIC PROC [in: Rope.ROPE, cellIndexStructure, inputPolynomialRing, inputSignatureStructure, minPolyRing, formulaAlgebra: AC.Structure] RETURNS [cells: QET.CellSeq] ~ {
stream: IO.STREAM ← IO.RIS[in];
termToken: Rope.ROPE;
[cells, termToken ] ← ReadCellSeq[stream, cellIndexStructure, inputPolynomialRing, inputSignatureStructure, minPolyRing, formulaAlgebra];
};
CellSeqToRope: PUBLIC PROC [in: QET.CellSeq] RETURNS [out: Rope.ROPE] ~ {
out ← "";
FOR i:NAT IN [1..in.lengthPlus1 - 1] DO
out ← Rope.Cat[out, CellToRope[in[i]],"\n" ];
ENDLOOP;
};
WriteCellSeq: PUBLIC PROC [in: QET.CellSeq, out: IO.STREAM] ~ {
cellSeqRope: Rope.ROPE ← CellSeqToRope[in];
out.PutF["\n %g \n", IO.rope[cellSeqRope] ];
};
Cell Indices and Cell Index Sequences
ReadCellIndex: PUBLIC PROC [in: IO.STREAM, cellIndexStructure: AC.Structure] RETURNS [out: QET.CellIndex] ~ {
RETURN[cellIndexStructure.class.read[in, cellIndexStructure] ];
};
CellIndexFromRope: PUBLIC PROC [in: Rope.ROPE, cellIndexStructure: AC.Structure] RETURNS [out: QET.CellIndex] = {
CIStream: IO.STREAM ← IO.RIS[in];
RETURN[ ReadCellIndex[ CIStream, cellIndexStructure ] ];
};
CellIndexToRope: PUBLIC PROC [in: QET.CellIndex] RETURNS [out: Rope.ROPE] ~ {
RETURN[in.structure.class.toRope[in] ];
};
WriteCellIndex: PUBLIC PROC [in: QET.CellIndex, out: IO.STREAM] = {
CIRope: Rope.ROPE ← CellIndexToRope[in];
out.PutF["%g\n", IO.rope[CIRope] ];
};
CellDimension: PUBLIC PROC [in: QET.CellIndex] RETURNS [dimension: CARDINAL] ~ {
inData: PTS.PointData ← NARROW[in.data];
dimension ← 0;
FOR i:NAT IN [1..inData.dimensionPlus1 - 1] DO
coodData: Ints.IntData ← NARROW[inData[i].data];
dimension ← dimension + coodData^ MOD 2;
ENDLOOP;
};
CellIndexEqual: PUBLIC PROC [in1, in2: QET.CellIndex] RETURNS [BOOL] = {
RETURN[in1.structure.class.equal[in1, in2] ];
};
ReadCellIndexSeq: PUBLIC PROC [in: IO.STREAM, cellIndexStructure: AC.Structure] RETURNS [out: QET.CellIndexSeq] ~ {
puncChar: CHAR;
nextIndex: QET.CellIndex;
length: NAT ← 0;
indexList, indexListTail: LIST OF QET.CellIndex ← NIL;
[]← in.SkipWhitespace[];
puncChar ← in.GetChar[];
IF puncChar # '( THEN ERROR;
WHILE puncChar # ') DO
nextIndex ← ReadCellIndex[in, cellIndexStructure];
length ← length + 1;
IF indexList=NIL THEN indexList ← indexListTail ←LIST[nextIndex] ELSE
{ indexListTail.rest ← LIST[nextIndex]; indexListTail ← indexListTail.rest };
[]← in.SkipWhitespace[];
puncChar ← in.GetChar[];
ENDLOOP;
out ← NEW[QET.CellIndexSeqRec[length]];
FOR i:NAT IN [1..length] DO
out[i] ← indexList.first;
indexList ← indexList.rest;
ENDLOOP;
};
CellIndexSeqFromRope: PUBLIC PROC [in: Rope.ROPE, cellIndexStructure: AC.Structure] RETURNS [out: QET.CellIndexSeq] ~ {
CISStream: IO.STREAM ← IO.RIS[in];
RETURN[ ReadCellIndexSeq[ CISStream, cellIndexStructure ] ];
};
CellIndexSeqToRope: PUBLIC PROC [in: QET.CellIndexSeq] RETURNS [out: Rope.ROPE] ~ {
out ← "(\n";
FOR i:NAT IN [1..in.lengthPlus1 - 1] DO
out ← Rope.Concat[ out, CellIndexToRope[in[i]]];
IF i < in.lengthPlus1 - 1 THEN
out ← Rope.Concat[ out, ",\n"]
ELSE
out ← Rope.Concat[ out, "\n"]
ENDLOOP;
out ← Rope.Concat[ out, ")\n" ];
};
WriteCellIndexSeq: PUBLIC PROC [in: QET.CellIndexSeq, out: IO.STREAM] ~ {
CISRope: Rope.ROPE ← CellIndexSeqToRope[in];
out.PutF["%g", IO.rope[CISRope] ];
};
Signatures
ReadSignature: PUBLIC PROC [in: IO.STREAM, signatureStructure: AC.Structure] RETURNS [out: QET.Signature] ~ {
RETURN[signatureStructure.class.read[in, signatureStructure] ];
};
SignatureToRope: PUBLIC PROC [in: QET.Signature] RETURNS [out: Rope.ROPE] ~ {
RETURN[in.structure.class.toRope[in] ];
};
Sample Points
ReadSamplePoint: PUBLIC PROC [in: IO.STREAM, inputPolynomialRing, minPolyRing: AC.Structure] RETURNS [out: QET.SamplePoint] ~ {
polyRingData: POL.PolynomialRingData ← NARROW[inputPolynomialRing.instanceData];
V: VARS.VariableSeq ← polyRingData.allVariables;
dimension: NAT ← V.lengthPlus1 - 1;
token: Rope.ROPE;
primitiveElement: AN.AlgebraicNumber;
basePoint: PTS.Point;
definingPolynomial: POL.Polynomial;
isolatingInterval: RI.RatInterval;
point: PTS.Point;
numberField, algebraicPolynomialsRing: AC.Structure;
token ← in.GetID[];
IF Rope.Equal[token, "NILSAMPLEPOINT"] THEN RETURN[NIL]
ELSE IF Rope.Equal[token, "EXTENDEDSAMPLEPOINT"] THEN {
basePointStructure: AC.Structure;
token ← in.GetID[];
IF NOT Rope.Equal[token, "PRIMITIVEELEMENT"] THEN ERROR;
primitiveElement ← AN.RealAlgebraicNumbers.class.read[in, AN.RealAlgebraicNumbers];
numberField ← EF.MakeExtensionField[primitiveElement];
token ← in.GetID[];
IF NOT Rope.Equal[token, "BASEPOINT"] THEN ERROR;
basePointStructure ← PTS.MakePointStructure[numberField, dimension - 1];
basePoint ← PTS.Read[in, basePointStructure];
token ← in.GetID[];
algebraicPolynomialsRing ← POL.MakePolynomialStructure[numberField, VARS.MainVariable[V] ];
IF NOT Rope.Equal[token, "DEFININGPOLYNOMIAL"] THEN ERROR;
definingPolynomial ← NARROW[algebraicPolynomialsRing.class.read[in, algebraicPolynomialsRing] ];
token ← in.GetID[];
IF NOT Rope.Equal[token, "ISOLATINGINTERVAL"] THEN ERROR;
isolatingInterval ← RI.RatIntervals.class.read[in, RI.RatIntervals];
out ← NEW[QET.SamplePointRec ← [
basePoint: basePoint,
definingPolynomial: definingPolynomial,
isolatingInterval: isolatingInterval
] ];
}
ELSE IF Rope.Equal[token, "PRIMITIVESAMPLEPOINT"] THEN {
pointStructure: AC.Structure;
token ← in.GetID[];
IF NOT Rope.Equal[token, "PRIMITIVEELEMENT"] THEN ERROR;
primitiveElement ← AN.RealAlgebraicNumbers.class.read[in, AN.RealAlgebraicNumbers];
numberField ← EF.MakeExtensionField[primitiveElement];
token ← in.GetID[];
IF NOT Rope.Equal[token, "POINT"] THEN ERROR;
pointStructure ← PTS.MakePointStructure[numberField, dimension];
point ← PTS.Read[in, pointStructure];
out ← NEW[QET.SamplePointRec ← [
point: point
] ];
}
ELSE ERROR;
};
SamplePointToRope: PUBLIC PROC [in: QET.SamplePoint] RETURNS [out: Rope.ROPE] ~ {
out ← "";
IF in=NIL THEN {
out ← Rope.Concat[out, "NILSAMPLEPOINT\n"];
RETURN;
};
IF in.basePoint # NIL THEN {
pointStructureData: PTS.PointStructureData ← NARROW[in.basePoint.structure.instanceData];
extensionFieldData: EF.ExtensionFieldData ← NARROW[pointStructureData.coordinateStructure.instanceData];
out ← Rope.Concat[out, "EXTENDEDSAMPLEPOINT\n"];
out ← Rope.Concat[out, "PRIMITIVEELEMENT\n"];
out ← Rope.Cat[out, AN.RealAlgebraicNumbers.class.toRope[extensionFieldData.primitiveElement], "\n"];
out ← Rope.Concat[out, "BASEPOINT\n"];
out ← Rope.Cat[out, PTS.ToRope[in.basePoint], "\n"];
out ← Rope.Concat[out, "DEFININGPOLYNOMIAL\n"];
out ← Rope.Cat[out, POL.PolyToRope[in.definingPolynomial], "\n"];
out ← Rope.Concat[out, "ISOLATINGINTERVAL\n"];
out ← Rope.Cat[out, RI.RatIntervals.class.toRope[in.isolatingInterval], "\n"];
}
ELSE {
pointStructureData: PTS.PointStructureData ← NARROW[in.point.structure.instanceData];
extensionFieldData: EF.ExtensionFieldData ← NARROW[pointStructureData.coordinateStructure.instanceData];
out ← Rope.Concat[out, "PRIMITIVESAMPLEPOINT\n"];
out ← Rope.Concat[out, "PRIMITIVEELEMENT\n"];
out ← Rope.Cat[out, AN.RealAlgebraicNumbers.class.toRope[extensionFieldData.primitiveElement], "\n"];
out ← Rope.Concat[out, "POINT\n"];
out ← Rope.Cat[out, PTS.ToRope[in.point], "\n"];
};
};
Adjacencies
ReadAdjacencies: PUBLIC PROC [in: IO.STREAM, cad: QET.Cad] ~ {
token: Rope.ROPE ← in.GetID[];
puncChar: CHAR;
dimension: CARDINAL;
inducedCad: QET.Cad;
cells: QET.CellSeq;
cellIndexStructure: AC.Structure;
IF Rope.Equal[token, "BEGINADJACENCIES"] THEN token ← in.GetID[]; -- flush if present
IF NOT Rope.Equal[token, "DIMENSION"] THEN ERROR;
dimension ← in.GetCard[];
inducedCad ← GetInducedCad[cad, dimension];
cellIndexStructure ← PTS.MakePointStructure[Ints.Ints, dimension];
cells ← inducedCad.cells;
[]← in.SkipWhitespace[];
puncChar ← in.GetChar[];
IF puncChar # '( THEN ERROR;
[]← in.SkipWhitespace[];
puncChar ← in.PeekChar[];
WHILE puncChar # ') DO
ReadAdjacency[in, cellIndexStructure, cells];
[]← in.SkipWhitespace[];
puncChar ← in.PeekChar[];
ENDLOOP;
puncChar ← in.GetChar[];
};
AdjacenciesFromRope: PUBLIC PROC [in: Rope.ROPE, cad: QET.Cad] ~ {
stream: IO.STREAM ← IO.RIS[in];
ReadAdjacencies[stream, cad];
};
ReadAdjacency: PUBLIC PROC [in: IO.STREAM, cellIndexStructure: AC.Structure, cells: QET.CellSeq] ~ {
puncChar: CHAR;
firstCell, secondCell: QET.Cell;
firstCellIndex, secondCellIndex: QET.CellIndex;
[]← in.SkipWhitespace[];
puncChar ← in.GetChar[ ];
IF puncChar # '( THEN ERROR;
firstCellIndex ← ReadCellIndex[in, cellIndexStructure];
[]← in.SkipWhitespace[];
puncChar ← in.GetChar[ ];
IF puncChar # ', THEN ERROR;
secondCellIndex ← ReadCellIndex[in, cellIndexStructure];
[]← in.SkipWhitespace[];
puncChar ← in.GetChar[ ];
IF puncChar # ') THEN ERROR;
firstCell ← LookupCell[cells, firstCellIndex];
secondCell ← LookupCell[cells, secondCellIndex];
AddAdjacency[firstCell, secondCell];
AddAdjacency[secondCell, firstCell];
};
AddAdjacency: PROC [first, second: QET.Cell] ~ {
firstAdjacencies: QET.CellSeq ← first.adjacentCells;
new: BOOL ← TRUE;
firstLengthPlus1: QET.CellSeqBoundsType;
IF firstAdjacencies = NIL THEN firstLengthPlus1 ← 1 ELSE {
firstLengthPlus1 ← firstAdjacencies.lengthPlus1;
FOR i:NAT IN [1..firstLengthPlus1-1] DO
IF firstAdjacencies[i] = second THEN { new ← FALSE; EXIT };
ENDLOOP;
};
IF new THEN {
newFirstAdjacencies: QET.CellSeq ← NEW[QET.CellSeqRec[firstLengthPlus1]];
FOR i:NAT IN [1..firstLengthPlus1-1] DO
newFirstAdjacencies[i] ← firstAdjacencies[i];
ENDLOOP;
newFirstAdjacencies[firstLengthPlus1] ← second;
first.adjacentCells ← newFirstAdjacencies;
};
};
CellAdjacenciesToRope: PUBLIC PROC [in: QET.Cell] RETURNS [out: Rope.ROPE] ~ {
out ← "";
IF in.adjacentCells = NIL THEN RETURN[out];
FOR i:NAT IN [1..in.adjacentCells.lengthPlus1 - 1] DO
out ← Rope.Cat[ out, "( " , CellIndexToRope[in.index], " , "];
out ← Rope.Cat[ out, CellIndexToRope[in.adjacentCells[i].index], " )\n"];
ENDLOOP;
};
CadAdjacenciesToRope: PUBLIC PROC [in: QET.Cad] RETURNS [out: Rope.ROPE] ~ {
cad: QET.Cad ← in;
out ← "";
FOR i:NAT DECREASING IN [1..in.dimension] DO
out ← Rope.Concat[out, "BEGINADJACENCIES\n"];
out ← Rope.Concat[out, "DIMENSION\n"];
out ← Rope.Cat[out, Convert.RopeFromCard[cad.dimension], "\n"];
out ← Rope.Concat[out, "(\n"];
IF cad.cells#NIL THEN FOR j:NAT IN [1..cad.cells.lengthPlus1 - 1] DO
out ← Rope.Cat[out, CellAdjacenciesToRope[cad.cells[j] ] ];
ENDLOOP;
out ← Rope.Concat[out, ")\n\n"];
cad ← cad.inducedCad;
ENDLOOP;
};
WriteCadAdjacencies: PUBLIC PROC [in: QET.Cad, out: IO.STREAM] ~ {
out.PutRope[ CadAdjacenciesToRope[in] ]
};
Covering Sets and Triangulations
ReadCoveringSets: PUBLIC PROC [in: IO.STREAM, cad: QET.Cad] ~ {
token: Rope.ROPE ← in.GetID[];
index: QET.CellIndex;
dimension: CARDINAL;
inducedCad: QET.Cad;
cellIndexStructure: AC.Structure;
cells: QET.CellSeq;
cell: QET.Cell;
IF Rope.Equal[token, "BEGINCOVERINGSETS"] THEN token ← in.GetID[]; -- flush if present
IF NOT Rope.Equal[token, "DIMENSION"] THEN ERROR;
dimension ← in.GetCard[];
inducedCad ← GetInducedCad[cad, dimension];
cellIndexStructure ← PTS.MakePointStructure[Ints.Ints, dimension];
cells ← inducedCad.cells;
token ← in.GetID[];
WHILE Rope.Equal[token, "CELLINDEX"] DO
index ← ReadCellIndex[in, cellIndexStructure];
cell ← LookupCell[cells, index];
token ← in.GetID[];
IF NOT Rope.Equal[token, "BEGINCOVERINGSET"] THEN ERROR;
cell.coveringSet ← ReadCellCoveringSet[in, dimension];
token ← in.GetID[];
ENDLOOP;
IF NOT Rope.Equal[token, "ENDCOVERINGSETS"] THEN ERROR;
};
CoveringSetsFromRope: PUBLIC PROC [in: Rope.ROPE, cad: QET.Cad] ~ {
stream: IO.STREAM ← IO.RIS[in];
ReadCoveringSets[stream, cad];
};
ReadCellCoveringSet: PUBLIC PROC [in: IO.STREAM, dimension: CARDINAL] RETURNS [out: QET.CoveringSet] ~ {
token: Rope.ROPE ← in.GetID[];
pointStructure: AC.Structure ← PTS.MakePointStructure[Reals.Reals, dimension];
pointList, pointListTail: LIST OF Points.Point ← NIL;
columnList, columnListTail: LIST OF QET.Column ← NIL;
column: QET.Column;
point: PTS.Point;
columnListLength: NAT ← 0;
WHILE Rope.Equal[token, "BEGINCOLUMN"] DO
pointListLength: NAT ← 0;
[]← in.SkipWhitespace[];
WHILE in.PeekChar[]='( DO
point ← PTS.Read[in, pointStructure];
pointListLength ← pointListLength + 1;
IF pointList=NIL THEN pointList ← pointListTail ←LIST[point] ELSE
pointListTail ← pointListTail.rest ← LIST[point];
[]← in.SkipWhitespace[];
ENDLOOP;
token ← in.GetID[];
IF NOT Rope.Equal[token, "ENDCOLUMN"] THEN ERROR;
column ← NEW[QET.ColumnRec[pointListLength]]; -- assumes length >=1
FOR i:NAT IN [1..pointListLength] DO
column[i] ← pointList.first;
pointList ← pointList.rest;
ENDLOOP;
columnListLength ← columnListLength + 1;
IF columnList=NIL THEN columnList ← columnListTail ←LIST[column] ELSE
columnListTail ← columnListTail.rest ← LIST[column];
token ← in.GetID[];
ENDLOOP;
IF NOT Rope.Equal[token, "ENDCOVERINGSET"] THEN ERROR;
out ← NEW[QET.CoveringSetRec[columnListLength]]; -- assumes length >=1
FOR i:NAT IN [1..columnListLength] DO
out.data[i] ← columnList.first;
columnList ← columnList.rest;
ENDLOOP;
RETURN[out];
};
CellCoveringSetToRope: PUBLIC PROC [in: QET.Cell] RETURNS [out: Rope.ROPE] ~ {
dimension: NAT ← in.dimension;
pointStructure: AC.Structure ← PTS.MakePointStructure[Reals.Reals, dimension];
coveringSet: QET.CoveringSet ← in.coveringSet;
IF coveringSet = NIL THEN RETURN[NIL];
out ← Rope.Concat[out, "CELLINDEX\n"];
out ← Rope.Cat[out, CellIndexToRope[in.index], "\n" ];
out ← Rope.Concat[out, "BEGINCOVERINGSET\n"];
FOR i:NAT IN [1..coveringSet.columnsPlusOne-1] DO
out ← Rope.Concat[out, "BEGINCOLUMN\n"];
FOR j:NAT IN [1..coveringSet[i].pointsPlusOne-1] DO
out ← Rope.Cat[out, PTS.ToRope[coveringSet[i][j] ], "\n"];
ENDLOOP;
out ← Rope.Concat[out, "ENDCOLUMN\n"];
ENDLOOP;
out ← Rope.Concat[out, "ENDCOVERINGSET\n"];
};
CadCoveringSetsToRope: PUBLIC PROC [in: QET.Cad] RETURNS [out: Rope.ROPE] ~ {
cad: QET.Cad ← in;
out ← "";
FOR i:NAT DECREASING IN [1..in.dimension] DO
out ← Rope.Concat[out, "BEGINCOVERINGSETS\n"];
out ← Rope.Concat[out, "DIMENSION\n"];
out ← Rope.Cat[out, Convert.RopeFromCard[cad.dimension], "\n"];
IF cad.cells#NIL THEN FOR j:NAT IN [1..cad.cells.lengthPlus1 - 1] DO
rope: Rope.ROPE ← CellCoveringSetToRope[cad.cells[j] ];
IF rope # NIL THEN out ← Rope.Cat[out, rope, "\n"];
ENDLOOP;
out ← Rope.Concat[out, "ENDCOVERINGSETS\n\n"];
cad ← cad.inducedCad;
ENDLOOP;
};
WriteCadCoveringSets: PUBLIC PROC [in: QET.Cad, out: IO.STREAM] ~ {
out.PutRope[ CadCoveringSetsToRope[in] ]
};
CoveringSetSize: PUBLIC PROC [in: QET.CoveringSet] RETURNS [size: CARDINAL] ~ {
size ← 0;
IF in = NIL THEN RETURN[size];
FOR i:NAT IN [1..in.columnsPlusOne-1] DO
size ← size + in[i].pointsPlusOne-1;
ENDLOOP;
};
CountTriangles: PUBLIC PROC [in: QET.CoveringSet] RETURNS [number: CARDINAL] ~ {
K, i, j: NAT;
M, I, J: CARDINAL ← 0;
IF in = NIL THEN RETURN[0];
IF in.columnsPlusOne-1 = 1 THEN RETURN[0];
number ← 0;
K ← 0;
WHILE K < in.columnsPlusOne-1 DO
K ← K + 1;
IF K = in.columnsPlusOne-1 THEN EXIT;
I ← in[K].pointsPlusOne-1;
J ← in[K+1].pointsPlusOne-1;
IF I < 2 AND J < 2 THEN {
M ← M + I;
LOOP;
};
IF I = 0 OR J = 0 THEN {
M ← M + I;
LOOP;
};
i ← 1; j ← 1;
WHILE i < I AND j < J DO
number ← number + 2;
i ← i+1;
j ← j+1;
ENDLOOP;
WHILE i < I DO
number ← number + 1;
i ← i+1;
ENDLOOP;
WHILE j < J DO
number ← number + 1;
j ← j+1;
ENDLOOP;
M ← M + I;
ENDLOOP;
};
GenerateTriangles: PUBLIC PROC [in: QET.CoveringSet] RETURNS [out: QET.TriangleSeq] ~ {
K, i, j: NAT;
M, I, J: CARDINAL ← 0;
numTriangles: CARDINAL ← CountTriangles[in];
number: CARDINAL ← 0;
IF numTriangles = 0 THEN RETURN[NIL];
out ← NEW[QET.TriangleSeqRec[numTriangles]];
K ← 0;
WHILE K < in.columnsPlusOne-1 DO
K ← K + 1;
IF K = in.columnsPlusOne-1 THEN EXIT;
I ← in[K].pointsPlusOne-1;
J ← in[K+1].pointsPlusOne-1;
IF I < 2 AND J < 2 THEN {
M ← M + I;
LOOP;
};
IF I = 0 OR J = 0 THEN {
M ← M + I;
LOOP;
};
i ← 1; j ← 1;
WHILE i < I AND j < J DO
number ← number + 1;
out[number].firstVertex ← M+i;
out[number].secondVertex ← M+i+1;
out[number].thirdVertex ← M+I+j;
i ← i+1;
number ← number + 1;
out[number].firstVertex ← M+I+j;
out[number].secondVertex ← M+i;
out[number].thirdVertex ← M+I+j+1;
j ← j+1;
ENDLOOP;
WHILE i < I DO
number ← number + 1;
out[number].firstVertex ← M+i;
out[number].secondVertex ← M+i+1;
out[number].thirdVertex ← M+I+j;
i ← i+1;
ENDLOOP;
WHILE j < J DO
number ← number + 1;
out[number].firstVertex ← M+I+j;
out[number].secondVertex ← M+i;
out[number].thirdVertex ← M+I+j+1;
j ← j+1;
ENDLOOP;
M ← M + I;
ENDLOOP;
};
CoveringSetToThreeDWorldRope: PUBLIC PROC [in: QET.CoveringSet] RETURNS [out: Rope.ROPE] ~ {
dimension: NAT;
out ← NIL;
IF in = NIL THEN RETURN[out];
IF in.columnsPlusOne >= 2 AND in[1].pointsPlusOne >= 2 THEN
dimension ← NARROW[in[1][1].data, PTS.PointData].dimensionPlus1 - 1;
FOR i:NAT IN [1..in.columnsPlusOne-1] DO
FOR j:NAT IN [1..in[i].pointsPlusOne-1] DO
FOR k:NAT IN [1..dimension] DO
cood: AC.Object ← NARROW[in[i][j].data, PTS.PointData][k];
value: REAL ← NARROW[cood.data, Reals.RealData]^;
out ← Rope.Cat[out, Convert.RopeFromReal[value], " "];
ENDLOOP;
out ← Rope.Concat[out,"\n"];
ENDLOOP;
ENDLOOP;
};
GenerateThreeDWorldTrianglesRope: PUBLIC PROC [in: QET.CoveringSet] RETURNS [out: Rope.ROPE] ~ {
K, i, j: NAT;
M, I, J: CARDINAL ← 0;
IF in = NIL THEN RETURN[NIL];
IF in.columnsPlusOne-1 = 1 THEN RETURN[NIL];
out ← NIL;
K ← 0;
WHILE K < in.columnsPlusOne-1 DO
K ← K + 1;
IF K = in.columnsPlusOne-1 THEN EXIT;
I ← in[K].pointsPlusOne-1;
J ← in[K+1].pointsPlusOne-1;
IF I < 2 AND J < 2 THEN {
M ← M + I;
LOOP;
};
IF I = 0 OR J = 0 THEN {
M ← M + I;
LOOP;
};
i ← 1; j ← 1;
WHILE i < I AND j < J DO
out ← Rope.Cat[out, "3 ", Convert.RopeFromCard[M+i], " "];
out ← Rope.Cat[out, Convert.RopeFromCard[M+i+1], " ", Convert.RopeFromCard[M+I+j], " \n"];
i ← i+1;
out ← Rope.Cat[out, "3 ", Convert.RopeFromCard[M+I+j], " "];
out ← Rope.Cat[out, Convert.RopeFromCard[M+i], " ", Convert.RopeFromCard[M+I+j+1], " \n"];
j ← j+1;
ENDLOOP;
WHILE i < I DO
out ← Rope.Cat[out, "3 ", Convert.RopeFromCard[M+i], " "];
out ← Rope.Cat[out, Convert.RopeFromCard[M+i+1], " ", Convert.RopeFromCard[M+I+j], " \n"];
i ← i+1;
ENDLOOP;
WHILE j < J DO
out ← Rope.Cat[out, "3 ", Convert.RopeFromCard[M+I+j], " "];
out ← Rope.Cat[out, Convert.RopeFromCard[M+i], " ", Convert.RopeFromCard[M+I+j+1], " \n"];
j ← j+1;
ENDLOOP;
M ← M + I;
ENDLOOP;
};
CellToThreeDWorldRope: PUBLIC PROC [in: QET.Cell] RETURNS [out: Rope.ROPE] ~ {
IF in.coveringSet = NIL THEN RETURN[NIL];
IF in.dimension # 2 THEN RETURN[NIL];
out ← Rope.Cat[Convert.RopeFromCard[CoveringSetSize[in.coveringSet]], " ", Convert.RopeFromCard[CountTriangles[in.coveringSet]], "\n"];
out ← Rope.Concat[out, CoveringSetToThreeDWorldRope[in.coveringSet] ];
out ← Rope.Cat[out, GenerateThreeDWorldTrianglesRope[in.coveringSet], "\n"];
};
CAD Files
ReadCadFile: PUBLIC PROC [in: IO.STREAM] RETURNS [out: QET.Cad] ~ {
token: Rope.ROPE;
out ← ReadCad[in];
[]← in.SkipWhitespace[];
IF NOT in.EndOf[] THEN token ← in.GetID[];
WHILE NOT in.EndOf[] DO
SELECT TRUE FROM
Rope.Equal[token, "BEGINADJACENCIES"] =>
ReadAdjacencies[in, out];
Rope.Equal[token, "BEGINCOVERINGSETS"] =>
ReadCoveringSets[in, out];
ENDCASE => ERROR;
[]← in.SkipWhitespace[];
IF NOT in.EndOf[] THEN token ← in.GetID[];
ENDLOOP;
RETURN[out];
};
CadFileFromRope: PUBLIC PROC [in: Rope.ROPE] RETURNS [out: QET.Cad] ~ {
stream: IO.STREAM ← IO.RIS[in];
out ← ReadCadFile[stream];
};
CadFileToRope: PUBLIC PROC [in: QET.Cad] RETURNS [out: Rope.ROPE] ~ {
out ← CadToRope[in];
out ← Rope.Concat[out, CadAdjacenciesToRope[in] ];
out ← Rope.Cat[out, CadCoveringSetsToRope[in] ];
};
WriteCadFile: PUBLIC PROC [in: QET.Cad, out: IO.STREAM] ~ {
cadRope: Rope.ROPE ← CadFileToRope[in];
out.PutRope[cadRope];
};
END.