Voici un proof of concept, écrit en OCaml (je n'ai pas de matlab sous la main).
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
| #!/usr/bin/ocaml
(* Vector type *)
type vec = float list
(* Vector operators *)
let vec_to_string = function
| x::y::z::t::w::[] -> Printf.sprintf "(%.2f, %.2f, %.2f, %.2f, %.2f)" x y z t w
| _ -> failwith "vec_to_string: bad format."
let vec_add =
List.map2 ( +. )
let vec_scal_prod r =
List.map ( ( *.) r)
let vec_outer_prod u v =
let w = List.map2 ( *. ) u v in
List.fold_left (+.) 0. w
let vec_norme u =
sqrt (vec_outer_prod u u)
(* Vector generator *)
let vec_nul () =
[0.;0.;0.;0.;0.]
let vec_random () =
List.map (fun _ -> Random.float 10.) (vec_nul ())
(* Algorithm functions *)
let compute_linear_combination a v =
let rec loop acc k =
if k < 0
then acc
else
let akvk = vec_scal_prod a.(k) v.(k) in
loop (vec_add acc akvk) (pred k)
in
loop (vec_nul ()) (Array.length a - 1)
let comute_error a v t =
let s = compute_linear_combination a v in
let minus_s = vec_scal_prod (-1.) s in
vec_add t minus_s
let compute_coefficient k a v t =
let a' = Array.mapi (fun i x -> if i = k then 0. else x) a in
let e = comute_error a' v t in
let sq x = x *. x in
let perfect_coef = (vec_outer_prod e v.(k)) /. (sq (vec_norme v.(k))) in
match perfect_coef with
| x when x < 0. -> 0.
| x when x > 1. -> 1.
| x -> x
let main () =
Random.self_init () ;
(* Make 42 random vectors *)
let v = Array.init 42 (fun _ -> vec_random ()) in
(* Initiate the 42 coef to zero (0) *)
let a = Array.make (Array.length v) 0. in
(* Make the target vector *)
let t = vec_random () in
(* Compute all the coefficients *)
let pass () =
let rec loop k =
if k >= (Array.length a)
then ()
else begin
a.(k) <- compute_coefficient k a v t ;
loop (succ k)
end
in
loop 0
in
(* make a new pass until target acquired or maxiter reached *)
let print_pass n =
Printf.printf "Pass #%d, Delta: %.3f\n" n (vec_norme (comute_error a v t))
in
let rec make_pass maxiter =
if (maxiter <= 0) || (vec_norme (comute_error a v t) < 1e-2)
then ()
else begin pass () ; print_pass maxiter ; make_pass (pred maxiter) end
in
make_pass 50 ;
print_pass 0
let _ = main () |
À l'exécution, en fonction des vecteurs générés au début du programme, on parvient à la précision souhaitée ou nous. Mais il est important de noter qu'à chaque passe, l'erreur diminue.
yscialom@Lucy:/tmp$ ./combi.ml
Pass #50, Delta: 7.092
Pass #49, Delta: 3.594
Pass #48, Delta: 2.331
Pass #47, Delta: 1.751
Pass #46, Delta: 1.364
Pass #45, Delta: 1.040
Pass #44, Delta: 0.783
Pass #43, Delta: 0.584
Pass #42, Delta: 0.431
Pass #41, Delta: 0.321
Pass #40, Delta: 0.251
Pass #39, Delta: 0.195
Pass #38, Delta: 0.150
Pass #37, Delta: 0.115
Pass #36, Delta: 0.086
Pass #35, Delta: 0.063
Pass #34, Delta: 0.045
Pass #33, Delta: 0.031
Pass #32, Delta: 0.022
Pass #31, Delta: 0.015
Pass #30, Delta: 0.010
Pass #29, Delta: 0.007
Pass #0, Delta: 0.007
Cdlt,
Chercher la combinaison linaire et satisfaction de contrainte
Répondre avec citation
Envoyé par SKone
Partager