1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
|
unit Unit1;
interface
uses
Windows, Messages, SysUtils, Variants, Classes, Graphics, Controls, Forms,
Dialogs, StdCtrls, Math;
type
TForm1 = class(TForm)
Button1: TButton;
Memo1: TMemo;
procedure Button1Click(Sender: TObject);
procedure FormClose(Sender: TObject; var Action: TCloseAction);
procedure FormActivate(Sender: TObject);
private
{ Déclarations privées }
public
{ Déclarations publiques }
end;
var
Form1: TForm1;
implementation
{$R *.dfm}
type _3 = array[1..3] of Double;
type _1ligne = record
xyz,
abc : _3;
end;
var T : TextFile;
St : String;
P : pointer;
function Extract ( S : string; var f : _1ligne ) : boolean;
var U : string;
FF : array[1..6] of Double absolute f;
status,
i,j : integer;
begin
Extract := false;
S := Trim(S) + ' ';
j:=0;
repeat
u:='';
i:=1;
while S[i] <> ' ' do
inc(i);
U := Trim(copy(S,1,i));
delete(S,1,i);
S:=Trim(S) + ' ';
if ( U = '' ) then
exit;
inc(j);
val( U, FF[j], status);
if ( status <> 0 ) then
exit; //
until ( j = 6);
Extract := true;
end;
var P_Sz : integer;
n : integer;
M33,
bj,
aij : array[1..3,1..3] of double;
procedure Make_mat_33;
var s : integer;
f : _1ligne;
p1 : pointer;
j,k : Byte;
begin
FillChar(aij,0,sizeof(aij));
p1 := P;
for s:=1 to n do
begin
f := _1ligne(p1^);
inc(longint(p1),sizeof(_1ligne));
for j:=1 to 3 do
for k:=1 to j do // k <= j
aij[j,k]:= aij[j,k]+ f.xyz[j]*f.xyz[k];
end;
aij[1,2]:=aij[2,1];
aij[1,3]:=aij[3,1];
aij[2,3]:=aij[3,2];
end;
procedure Make_XK;
var s : integer;
f : _1ligne;
p1 : pointer;
j,b : Byte;
begin
FillChar(bj,0,sizeof(bj));
p1 := P;
for s:=1 to n do
begin
f := _1ligne(p1^);
inc(longint(p1),sizeof(_1ligne));
for b:=1 to 3 do
for j:=1 to 3 do
bj[j,b]:= bj[j,b]+ f.xyz[j]*f.abc[b];
end;
end;
function deter ( u11,u21,u31, u12,u22,u32, u13,u23,u33 : double ) : double;
begin
deter := u11 * (u22*u33-u32*u23) - u12*(u21*u33-u23*u31) + u13 * (u21*u32-u22*u31);
end;
var D, Residu_Max : double;
var Dx, Dy, Dz, R : array[1..3] of double;
const Error_Max_allowed = 1e-3;
function Control_des_ecarts : boolean;
var i : integer;
p1 : pointer;
f : _1ligne;
b : Byte;
begin
//
p1 := P;
Residu_Max:=0;
for i:=1 to n do
begin
f := _1ligne(p1^);
for b:=1 to 3 do
begin
R[b]:= f.abc[b] - M33[b,1]*f.xyz[1] - M33[b,2]*f.xyz[2] - M33[b,3] * f.xyz[3];
Residu_Max := Math.Max(Residu_Max,R[b]);
end;
inc(LongInt(p1),Sizeof(f));
end;
Control_des_ecarts := ( Residu_Max < Error_Max_allowed );
end;
function Make_Mij : boolean;
var b : byte; RR : string;
begin
(*
ici ki = a,b,c suivant i ki = Mi1*x + Mi2*y + Mi3*z à résoudre
en fait numériquement cela n'est jamais le cas
chaque equation 1..n presente un ecatr Epsilon[j] = [ki - (Mi1*x + Mi2*y + Mi3*z)] j= 1..n
on cherche à minimiser le Résidu sigma ( Espilon[j] ^2 )
R = Sigma [ ki(j) - Mi1* X(j) + Mi2 * Y(j) + Mi2 * Z(j) ] ^2
j = 1..n
soit @R/@Mia = 0 a=1..3
@R@Mia = Sigma [ ki(j) - Mi1* X(j) + Mi2 * Y(j) + Mi2 * Z(j) ] * 2 * HH[j] HH = X, Y ou Z pur a = 1,2 ,3
j = 1..n
donc le système résolvant s'écrit
Mi1 * Sigm X^2 + Mi2 * Sima ( X*Y ) + Mi3 * Sigma(X*Z) = Sima( X* Ki)
Mi1 * Sigm Y*X + Mi2 * Sima ( Y^2 ) + Mi3 * Sigma(Y*Z) = Sima( Y* Ki)
Mi1 * Sigm Z*X + Mi2 * Sima ( Z*Y ) + Mi3 * Sigma(Z^2) = Sima( Z* Ki)
on note que la matrice carré résolvant le système est la même pour i=1,2,3
=A> on la cacule qu'une X.
les solutions Mia diffèrent pour les differents i que par Sigma(x,y,z * Ki (a,b,c suivant i )
de + A est symétrique
*)
Make_mat_33;
D := deter( aij[1,1], aij[2,1],aij[3,1],
aij[1,2], aij[2,2],aij[3,2],
aij[1,3], aij[2,3],aij[3,3]
);
if ( abs( D ) < 1e-15 ) then
begin
Make_Mij := false;
ShowMessage('M n''existe pas système dégénéré...');
exit;
end;
Make_Mij := true;
D := 1/D;
Make_XK ;
for b:=1 to 3 do
begin
Dx[b] := deter( bj[1,b], aij[2,1],aij[3,1],
bj[2,b], aij[2,2],aij[3,2],
bj[3,b], aij[2,3],aij[3,3]
);
Dy[b] := deter( aij[1,1], bj[1,b],aij[3,1],
aij[1,2], bj[2,b],aij[3,2],
aij[1,3], bj[3,b],aij[3,3]
);
Dz[b] := deter( aij[1,1], aij[2,1],bj[1,b],
aij[1,2], aij[2,2],bj[2,b],
aij[1,3], aij[2,3],bj[3,b]
);
M33[b,1]:=D*Dx[b];
M33[b,2]:=D*Dy[b];
M33[b,3]:=D*Dz[b];
end;
str(Residu_Max:20:10,RR);
RR :=Trim(RR);
if ( Control_des_ecarts ) then
Form1.Caption := 'Matrice trouvée avec R max = ' + RR
else
ShowMessage('Aucune matrice 3X3 trouvée dans le cadre de la tolérance envisagée Résidu = ' + RR);
end;
function Make_M3x3 : boolean;
var i : byte;
begin
(*
pour chaque ligne on a ( k1,k2,k3 ) = ( a,b,c )
k1 M11 M12 M13 x
k2 = M21 M22 M23 * y
k3 M31 M32 M33 z
soit
ki = Mi1*x + Mi2*y + Mi3*z i=1..3
*)
if ( not Make_Mij ) then
begin
Make_M3x3 := false;
exit;
end;
Make_M3x3 := true;
end;
procedure TForm1.Button1Click(Sender: TObject);
var f : _1ligne;
fl : file;
i,j : byte;
coef : array[1..3,1..3] of String;
begin
Memo1.Clear;
if ( P <> NIL) then
begin
freemem(P,P_Sz);
P:=NIL;
end;
if ( FileExists('data.txt') ) then
begin
n := 0;
assignfile(fl,'data.bin');
rewrite(fl,sizeof(_1ligne));
assignfile(T,'data.txt');
reset(T);
// recherche des lignes acceptées => suite de 6 float ecrits sous forme ---.--- séparés par space
while not eof(T) do
begin
readln(T,St);
if ( Extract(St,f) ) then
begin
blockwrite(fl,f,1);
inc(n); // 1 de + OK !!
end;
end;
closefile(T);
closefile(fl);
if ( n < 3 ) then
begin
P := NIL;
erase(fl);
ShowMessage('Un minimum de 3 lignes valides sont NECESSAIRES !!!')
end
else
begin
P_Sz := n*Sizeof(_1ligne);
getmem(P,P_Sz);
reset(fl,P_Sz);
blockread(fl,p^,1);
closefile(fl);
erase(fl);
if ( Make_M3x3 ) then with Memo1 do
begin
for i:=1 to 3 do for j:=1 to 3 do
str(M33[i,j]:20:6,coef[i,j]);
lines.add ( Coef[1,1] + ' ' + Coef[1,2] + ' ' + Coef[1,3]);
lines.add ( Coef[2,1] + ' ' + Coef[2,2] + ' ' + Coef[2,3]);
lines.add ( Coef[3,1] + ' ' + Coef[3,2] + ' ' + Coef[3,3]);
end;
end
end
else
ShowMessage(' Le fichier text data.txt est absent');
end;
procedure TForm1.FormClose(Sender: TObject; var Action: TCloseAction);
begin
if ( P <> NIL) then
begin
freemem(P,P_Sz);
P:=NIL;
end;
end;
procedure TForm1.FormActivate(Sender: TObject);
begin
Memo1.Clear;
end;
initialization
P_Sz := 0;
P := NIL;
end. |
