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 | with Ada.Numerics.Elementary_Functions;
package body Prime_And_Coprime_Numbers with
Refined_State => (Prime_Data => Set.Is_Prime)
is
No_Number_Error : exception;
-----------------------------------------------------------------------------
package Set is
Is_Prime : Number_List_Type;
end Set;
function Valid_Prime_Data return Boolean is
(for all V in Number_List_Type'Range => Set.Is_Prime (V) = Is_Prime (V));
function Has_Prime (Low, High : Value_Type) return Boolean is
(for some V in Low .. High => Is_Prime (V));
-----------------------------------------------------------------------------
function Initialize_Coprime_List
(Value : in Value_Type)
return Number_List_Type
is
Result : Number_List_Type;
function Euclid
(Left : in Value_Type;
Right : in Value_Type)
return Value_Type
with
Pre => Left >= 2 and Right >= 0,
Contract_Cases =>
(Are_Coprime (Left, Right) => Euclid'Result = 1,
others => Euclid'Result > 1)
is
A : Value_Type range 1 .. Value_Type'Last;
B : Value_Type range 1 .. Value_Type'Last;
R : Value_Type range 0 .. Value_Type'Last;
begin
if Left = 0
then
return Right;
elsif Right = 0
then
return Left;
end if;
A := Left;
B := Right;
loop
pragma Loop_Invariant (A > 0 and B > 0);
pragma Loop_Invariant (not (A = 1 and B = 1));
pragma Loop_Invariant
(Are_Coprime (A, B) = Are_Coprime (Left, Right));
R := A mod B;
exit when R = 0;
pragma Assume (Are_Coprime (A, B) = Are_Coprime (B, R));
A := B;
B := R;
end loop;
pragma Assert (not (A = 1 and B > 0));
return B;
end Euclid;
begin
for Index in Result'Range
loop
Result (Index) := Euclid (Value, Index) = 1;
pragma Loop_Invariant
(for all V in Result'First .. Index =>
Result (V) = Are_Coprime (Value, V));
pragma Annotate (GNATprove, False_Positive, """Result"" might not be initialized",
"Result accessed in loop invariant only for indexes previously initialized");
end loop;
return Result;
end Initialize_Coprime_List;
-----------------------------------------------------------------------------
function Nearest_Number
(Number_List : in Number_List_Type;
Mode : in Nearest_Mode;
Value : in Value_Type)
return Value_Type
is
Right : Value_Type'Base := Value_Type'First;
Left : Value_Type'Base := Value_Type'First;
Right_Is_Out : Boolean;
Left_Is_Out : Boolean;
begin
if Number_List (Value)
then
return Value;
end if;
-- Search a number larger or equal to Value
if Mode = Up or else Mode = Absolute
then
Right := Value_Type'Succ (Value);
loop
Right_Is_Out := Right > Number_List'Last;
exit when
Right_Is_Out
or else Number_List (Right);
pragma Loop_Invariant (Right in Value + 1 .. Number_List'Last);
pragma Loop_Invariant
(for all V in Value .. Right => Number_List (V) = False);
Right := Value_Type'Succ (Right);
end loop;
else
Right_Is_Out := True;
end if;
-- Search a number lower than Value
if Mode = Down or else Mode = Absolute
then
Left := Value_Type'Pred (Value);
loop
Left_Is_Out := Left < Number_List'First;
exit when
Left_Is_Out
or else Number_List (Left);
pragma Loop_Invariant (Left in Number_List'First .. Value - 1);
pragma Loop_Invariant
(for all V in Left .. Value => Number_List (V) = False);
Left := Value_Type'Pred (Left);
end loop;
else
Left_Is_Out := True;
end if;
-- Select nearest prime
if not Right_Is_Out
then
if not Left_Is_Out
then
if abs (Left - Value) <= abs (Right - Value)
then
return Left;
else
return Right;
end if;
else
return Right;
end if;
else
if not Left_Is_Out
then
return Left;
else
raise No_Number_Error;
end if;
end if;
end Nearest_Number;
-----------------------------------------------------------------------------
function Nearest_Prime_Number
(Value : in Value_Type;
Mode : in Nearest_Mode)
return Value_Type
is
begin
return Nearest_Number (Set.Is_Prime, Mode, Value);
end Nearest_Prime_Number;
-----------------------------------------------------------------------------
procedure Eratosthenes with
Global => (Output => Set.Is_Prime),
Post => Valid_Prime_Data
is
Index_1 : Value_Type;
Index_3 : Value_Type'Base;
use Ada.Numerics.Elementary_Functions;
begin
Set.Is_Prime := (others => True);
Set.Is_Prime (0) := False;
Set.Is_Prime (1) := False;
for Index_2 in 2 .. Value_Type (Sqrt (Float (Max_Value)))
loop
pragma Loop_Invariant (Index_2 in 2 .. Max_Value);
pragma Loop_Invariant
(for all V in 0 .. Index_2 => Set.Is_Prime (V) = Is_Prime (V));
pragma Loop_Invariant
(for all V in Index_2 .. Set.Is_Prime'Last =>
Set.Is_Prime (V) =
(for all Div in 2 .. Index_2 - 1 => V mod Div /= 0));
if Set.Is_Prime (Index_2)
then
Index_1 := Index_2;
Index_3 := 2 * Index_1;
while Index_3 <= Max_Value
loop
Set.Is_Prime (Index_3) := False;
pragma Loop_Invariant
(for all V in 0 .. Index_2 => Set.Is_Prime (V) = Is_Prime (V));
pragma Loop_Invariant
(for all V in Index_2 .. Set.Is_Prime'Last =>
Set.Is_Prime (V) =
((for all Div in 2 .. Index_2 - 1 => V mod Div /= 0)
and
(if V in Index_2 + 1 .. Index_3 then V mod Index_2 /= 0)));
pragma Loop_Invariant (Index_3 in Index_2 .. Max_Value);
pragma Loop_Invariant (Index_3 mod Index_2 = 0);
Index_3 := Index_3 + Index_1;
end loop;
end if;
end loop;
end Eratosthenes;
-----------------------------------------------------------------------------
begin
Eratosthenes;
end Prime_And_Coprime_Numbers;
|