![[BACK]](/icons/cvsweb/back.gif){kind=icon}
Return to projective_transformations.adb CVS log ![[TXT]](/icons/cvsweb/text.gif){kind=icon}
![[DIR]](/icons/cvsweb/dir.gif){kind=icon}
Up to [local] / OpenXM_contrib / PHC / Ada / Homotopy|
Return to projective_transformations.adb CVS log | Up to [local] / OpenXM_contrib / PHC / Ada / Homotopy |
1.1 maekawa 1: with Standard_Floating_Numbers; use Standard_Floating_Numbers;
2: with Standard_Complex_Numbers; use Standard_Complex_Numbers;
3: with Standard_Natural_Vectors;
4:
5: package body Projective_Transformations is
6:
7: function Projective_Transformation ( p : Poly ) return Poly is
8:
9: deg : integer := Degree(p);
10: htd : Degrees
11: := new Standard_Natural_Vectors.Vector(1..Number_of_Unknowns(p)+1);
12: res : Poly := Null_Poly;
13:
14: procedure Embed_Term ( t : in Term; continue : out boolean ) is
15:
16: ht : Term;
17: sum : natural := 0;
18:
19: begin
20: ht.cf := t.cf;
21: for i in t.dg'range loop
22: sum := sum + t.dg(i);
23: htd(i) := t.dg(i);
24: end loop;
25: htd(htd'last) := deg-sum;
26: ht.dg := htd;
27: Add(res,ht);
28: continue := true;
29: end Embed_Term;
30: procedure Embed_Terms is new Visiting_Iterator(Embed_Term);
31:
32: begin
33: Embed_Terms(p);
34: Clear(htd);
35: return res;
36: end Projective_Transformation;
37:
38: procedure Projective_Transformation ( p : in out Poly ) is
39:
40: res : Poly := Projective_Transformation(p);
41:
42: begin
43: Copy(res,p); Clear(res);
44: end Projective_Transformation;
45:
46: function Projective_Transformation ( p : Poly_Sys ) return Poly_Sys is
47:
48: res : Poly_Sys(p'range);
49:
50: begin
51: for k in p'range loop
52: res(k) := Projective_Transformation(p(k));
53: end loop;
54: return res;
55: end Projective_Transformation;
56:
57: procedure Projective_Transformation ( p : in out Poly_Sys ) is
58: begin
59: for k in p'range loop
60: Projective_Transformation(p(k));
61: end loop;
62: end Projective_Transformation;
63:
64: function Projective_Transformation ( s : Solution ) return Solution is
65:
66: n : natural := s.n;
67: r : Solution(n+1);
68:
69: begin
70: r.v(1..n) := s.v(1..n);
71: r.v(n+1) := Create(1.0);
72: r.t := s.t;
73: r.m := s.m;
74: r.err := s.err;
75: r.rco := s.rco;
76: r.res := s.res;
77: return r;
78: end Projective_Transformation;
79:
80: function Projective_Transformation
81: ( sols : Solution_List ) return Solution_List is
82:
83: res,res_last : Solution_List;
84: tmp : Solution_List := sols;
85: ls : Link_to_Solution;
86:
87: begin
88: while not Is_Null(tmp) loop
89: ls := Head_Of(tmp);
90: Append(res,res_last,Projective_Transformation(ls.all));
91: tmp := Tail_Of(tmp);
92: end loop;
93: return res;
94: end Projective_Transformation;
95:
96: procedure Projective_Transformation ( sols : in out Solution_List ) is
97: begin
98: if Is_Null(sols)
99: then null;
100: else declare
101: temp : Solution_List := sols;
102: n : natural := Head_Of(sols).n;
103: l : Link_To_Solution;
104: s : Solution(n);
105: s2 : Solution(n+1);
106: begin
107: while not Is_Null(temp) loop
108: l := Head_Of(temp);
109: s := l.all;
110: s2 := Projective_Transformation(s);
111: Clear(l);
112: l := new Solution'(s2);
113: Set_Head(temp,l);
114: temp := Tail_Of(temp);
115: end loop;
116: end;
117: end if;
118: end Projective_Transformation;
119:
120: function Affine_Transformation ( s : Solution ) return Solution is
121:
122: n : natural := s.n;
123: r : Solution(n-1);
124: absvn : double_float := AbsVal(s.v(n));
125:
126: begin
127: for i in 1..(n-1) loop
128: if absvn + Create(1.0) = Create(1.0)
129: then r.v(i) := Create(10.0**10);
130: else r.v(i) := s.v(i) / s.v(n);
131: end if;
132: end loop;
133: r.t := s.t;
134: r.m := s.m;
135: r.err := s.err;
136: r.rco := s.rco;
137: r.res := s.res;
138: return r;
139: exception
140: when numeric_error => r.v(1..(n-1)) := (1..(n-1) => Create(10.0**10));
141: return r;
142: end Affine_Transformation;
143:
144: procedure Affine_Transformation ( sols : in out Solution_List ) is
145: begin
146: if Is_Null(sols)
147: then null;
148: else declare
149: n : natural := Head_Of(sols).n;
150: s1 : Solution(n);
151: s2 : Solution(n-1);
152: temp : Solution_List := sols;
153: l : Link_To_Solution;
154: begin
155: while not Is_Null(temp) loop
156: l := Head_Of(temp);
157: s1 := l.all;
158: s2 := Affine_Transformation(s1);
159: Clear(l);
160: l := new Solution'(s2);
161: Set_Head(temp,l);
162: temp := Tail_Of(temp);
163: end loop;
164: end;
165: end if;
166: end Affine_Transformation;
167:
168: end Projective_Transformations;