Coverage for /usr/lib/python3/dist-packages/scipy/optimize/_shgo_lib/_vertex.py: 22%

237 statements  

« prev     ^ index     » next       coverage.py v7.9.1, created at 2025-06-14 15:55 +0200

1import collections 

2from abc import ABC, abstractmethod 

3 

4import numpy as np 

5 

6from scipy._lib._util import MapWrapper 

7 

8 

9class VertexBase(ABC): 

10 """ 

11 Base class for a vertex. 

12 """ 

13 def __init__(self, x, nn=None, index=None): 

14 """ 

15 Initiation of a vertex object. 

16 

17 Parameters 

18 ---------- 

19 x : tuple or vector 

20 The geometric location (domain). 

21 nn : list, optional 

22 Nearest neighbour list. 

23 index : int, optional 

24 Index of vertex. 

25 """ 

26 self.x = x 

27 self.hash = hash(self.x) # Save precomputed hash 

28 

29 if nn is not None: 

30 self.nn = set(nn) # can use .indexupdate to add a new list 

31 else: 

32 self.nn = set() 

33 

34 self.index = index 

35 

36 def __hash__(self): 

37 return self.hash 

38 

39 def __getattr__(self, item): 

40 if item not in ['x_a']: 

41 raise AttributeError(f"{type(self)} object has no attribute " 

42 f"'{item}'") 

43 if item == 'x_a': 

44 self.x_a = np.array(self.x) 

45 return self.x_a 

46 

47 @abstractmethod 

48 def connect(self, v): 

49 raise NotImplementedError("This method is only implemented with an " 

50 "associated child of the base class.") 

51 

52 @abstractmethod 

53 def disconnect(self, v): 

54 raise NotImplementedError("This method is only implemented with an " 

55 "associated child of the base class.") 

56 

57 def star(self): 

58 """Returns the star domain ``st(v)`` of the vertex. 

59 

60 Parameters 

61 ---------- 

62 v : 

63 The vertex ``v`` in ``st(v)`` 

64 

65 Returns 

66 ------- 

67 st : set 

68 A set containing all the vertices in ``st(v)`` 

69 """ 

70 self.st = self.nn 

71 self.st.add(self) 

72 return self.st 

73 

74 

75class VertexScalarField(VertexBase): 

76 """ 

77 Add homology properties of a scalar field f: R^n --> R associated with 

78 the geometry built from the VertexBase class 

79 """ 

80 

81 def __init__(self, x, field=None, nn=None, index=None, field_args=(), 

82 g_cons=None, g_cons_args=()): 

83 """ 

84 Parameters 

85 ---------- 

86 x : tuple, 

87 vector of vertex coordinates 

88 field : callable, optional 

89 a scalar field f: R^n --> R associated with the geometry 

90 nn : list, optional 

91 list of nearest neighbours 

92 index : int, optional 

93 index of the vertex 

94 field_args : tuple, optional 

95 additional arguments to be passed to field 

96 g_cons : callable, optional 

97 constraints on the vertex 

98 g_cons_args : tuple, optional 

99 additional arguments to be passed to g_cons 

100 

101 """ 

102 super().__init__(x, nn=nn, index=index) 

103 

104 # Note Vertex is only initiated once for all x so only 

105 # evaluated once 

106 # self.feasible = None 

107 

108 # self.f is externally defined by the cache to allow parallel 

109 # processing 

110 # None type that will break arithmetic operations unless defined 

111 # self.f = None 

112 

113 self.check_min = True 

114 self.check_max = True 

115 

116 def connect(self, v): 

117 """Connects self to another vertex object v. 

118 

119 Parameters 

120 ---------- 

121 v : VertexBase or VertexScalarField object 

122 """ 

123 if v is not self and v not in self.nn: 

124 self.nn.add(v) 

125 v.nn.add(self) 

126 

127 # Flags for checking homology properties: 

128 self.check_min = True 

129 self.check_max = True 

130 v.check_min = True 

131 v.check_max = True 

132 

133 def disconnect(self, v): 

134 if v in self.nn: 

135 self.nn.remove(v) 

136 v.nn.remove(self) 

137 

138 # Flags for checking homology properties: 

139 self.check_min = True 

140 self.check_max = True 

141 v.check_min = True 

142 v.check_max = True 

143 

144 def minimiser(self): 

145 """Check whether this vertex is strictly less than all its 

146 neighbours""" 

147 if self.check_min: 

148 self._min = all(self.f < v.f for v in self.nn) 

149 self.check_min = False 

150 

151 return self._min 

152 

153 def maximiser(self): 

154 """ 

155 Check whether this vertex is strictly greater than all its 

156 neighbours. 

157 """ 

158 if self.check_max: 

159 self._max = all(self.f > v.f for v in self.nn) 

160 self.check_max = False 

161 

162 return self._max 

163 

164 

165class VertexVectorField(VertexBase): 

166 """ 

167 Add homology properties of a scalar field f: R^n --> R^m associated with 

168 the geometry built from the VertexBase class. 

169 """ 

170 

171 def __init__(self, x, sfield=None, vfield=None, field_args=(), 

172 vfield_args=(), g_cons=None, 

173 g_cons_args=(), nn=None, index=None): 

174 super().__init__(x, nn=nn, index=index) 

175 

176 raise NotImplementedError("This class is still a work in progress") 

177 

178 

179class VertexCacheBase: 

180 """Base class for a vertex cache for a simplicial complex.""" 

181 def __init__(self): 

182 

183 self.cache = collections.OrderedDict() 

184 self.nfev = 0 # Feasible points 

185 self.index = -1 

186 

187 def __iter__(self): 

188 for v in self.cache: 

189 yield self.cache[v] 

190 return 

191 

192 def size(self): 

193 """Returns the size of the vertex cache.""" 

194 return self.index + 1 

195 

196 def print_out(self): 

197 headlen = len(f"Vertex cache of size: {len(self.cache)}:") 

198 print('=' * headlen) 

199 print(f"Vertex cache of size: {len(self.cache)}:") 

200 print('=' * headlen) 

201 for v in self.cache: 

202 self.cache[v].print_out() 

203 

204 

205class VertexCube(VertexBase): 

206 """Vertex class to be used for a pure simplicial complex with no associated 

207 differential geometry (single level domain that exists in R^n)""" 

208 def __init__(self, x, nn=None, index=None): 

209 super().__init__(x, nn=nn, index=index) 

210 

211 def connect(self, v): 

212 if v is not self and v not in self.nn: 

213 self.nn.add(v) 

214 v.nn.add(self) 

215 

216 def disconnect(self, v): 

217 if v in self.nn: 

218 self.nn.remove(v) 

219 v.nn.remove(self) 

220 

221 

222class VertexCacheIndex(VertexCacheBase): 

223 def __init__(self): 

224 """ 

225 Class for a vertex cache for a simplicial complex without an associated 

226 field. Useful only for building and visualising a domain complex. 

227 

228 Parameters 

229 ---------- 

230 """ 

231 super().__init__() 

232 self.Vertex = VertexCube 

233 

234 def __getitem__(self, x, nn=None): 

235 try: 

236 return self.cache[x] 

237 except KeyError: 

238 self.index += 1 

239 xval = self.Vertex(x, index=self.index) 

240 # logging.info("New generated vertex at x = {}".format(x)) 

241 # NOTE: Surprisingly high performance increase if logging 

242 # is commented out 

243 self.cache[x] = xval 

244 return self.cache[x] 

245 

246 

247class VertexCacheField(VertexCacheBase): 

248 def __init__(self, field=None, field_args=(), g_cons=None, g_cons_args=(), 

249 workers=1): 

250 """ 

251 Class for a vertex cache for a simplicial complex with an associated 

252 field. 

253 

254 Parameters 

255 ---------- 

256 field : callable 

257 Scalar or vector field callable. 

258 field_args : tuple, optional 

259 Any additional fixed parameters needed to completely specify the 

260 field function 

261 g_cons : dict or sequence of dict, optional 

262 Constraints definition. 

263 Function(s) ``R**n`` in the form:: 

264 g_cons_args : tuple, optional 

265 Any additional fixed parameters needed to completely specify the 

266 constraint functions 

267 workers : int optional 

268 Uses `multiprocessing.Pool <multiprocessing>`) to compute the field 

269 functions in parrallel. 

270 

271 """ 

272 super().__init__() 

273 self.index = -1 

274 self.Vertex = VertexScalarField 

275 self.field = field 

276 self.field_args = field_args 

277 self.wfield = FieldWrapper(field, field_args) # if workers is not 1 

278 

279 self.g_cons = g_cons 

280 self.g_cons_args = g_cons_args 

281 self.wgcons = ConstraintWrapper(g_cons, g_cons_args) 

282 self.gpool = set() # A set of tuples to process for feasibility 

283 

284 # Field processing objects 

285 self.fpool = set() # A set of tuples to process for scalar function 

286 self.sfc_lock = False # True if self.fpool is non-Empty 

287 

288 self.workers = workers 

289 self._mapwrapper = MapWrapper(workers) 

290 

291 if workers == 1: 

292 self.process_gpool = self.proc_gpool 

293 if g_cons is None: 

294 self.process_fpool = self.proc_fpool_nog 

295 else: 

296 self.process_fpool = self.proc_fpool_g 

297 else: 

298 self.process_gpool = self.pproc_gpool 

299 if g_cons is None: 

300 self.process_fpool = self.pproc_fpool_nog 

301 else: 

302 self.process_fpool = self.pproc_fpool_g 

303 

304 def __getitem__(self, x, nn=None): 

305 try: 

306 return self.cache[x] 

307 except KeyError: 

308 self.index += 1 

309 xval = self.Vertex(x, field=self.field, nn=nn, index=self.index, 

310 field_args=self.field_args, 

311 g_cons=self.g_cons, 

312 g_cons_args=self.g_cons_args) 

313 

314 self.cache[x] = xval # Define in cache 

315 self.gpool.add(xval) # Add to pool for processing feasibility 

316 self.fpool.add(xval) # Add to pool for processing field values 

317 return self.cache[x] 

318 

319 def __getstate__(self): 

320 self_dict = self.__dict__.copy() 

321 del self_dict['pool'] 

322 return self_dict 

323 

324 def process_pools(self): 

325 if self.g_cons is not None: 

326 self.process_gpool() 

327 self.process_fpool() 

328 self.proc_minimisers() 

329 

330 def feasibility_check(self, v): 

331 v.feasible = True 

332 for g, args in zip(self.g_cons, self.g_cons_args): 

333 # constraint may return more than 1 value. 

334 if np.any(g(v.x_a, *args) < 0.0): 

335 v.f = np.inf 

336 v.feasible = False 

337 break 

338 

339 def compute_sfield(self, v): 

340 """Compute the scalar field values of a vertex object `v`. 

341 

342 Parameters 

343 ---------- 

344 v : VertexBase or VertexScalarField object 

345 """ 

346 try: 

347 v.f = self.field(v.x_a, *self.field_args) 

348 self.nfev += 1 

349 except AttributeError: 

350 v.f = np.inf 

351 # logging.warning(f"Field function not found at x = {self.x_a}") 

352 if np.isnan(v.f): 

353 v.f = np.inf 

354 

355 def proc_gpool(self): 

356 """Process all constraints.""" 

357 if self.g_cons is not None: 

358 for v in self.gpool: 

359 self.feasibility_check(v) 

360 # Clean the pool 

361 self.gpool = set() 

362 

363 def pproc_gpool(self): 

364 """Process all constraints in parallel.""" 

365 gpool_l = [] 

366 for v in self.gpool: 

367 gpool_l.append(v.x_a) 

368 

369 G = self._mapwrapper(self.wgcons.gcons, gpool_l) 

370 for v, g in zip(self.gpool, G): 

371 v.feasible = g # set vertex object attribute v.feasible = g (bool) 

372 

373 def proc_fpool_g(self): 

374 """Process all field functions with constraints supplied.""" 

375 for v in self.fpool: 

376 if v.feasible: 

377 self.compute_sfield(v) 

378 # Clean the pool 

379 self.fpool = set() 

380 

381 def proc_fpool_nog(self): 

382 """Process all field functions with no constraints supplied.""" 

383 for v in self.fpool: 

384 self.compute_sfield(v) 

385 # Clean the pool 

386 self.fpool = set() 

387 

388 def pproc_fpool_g(self): 

389 """ 

390 Process all field functions with constraints supplied in parallel. 

391 """ 

392 self.wfield.func 

393 fpool_l = [] 

394 for v in self.fpool: 

395 if v.feasible: 

396 fpool_l.append(v.x_a) 

397 else: 

398 v.f = np.inf 

399 F = self._mapwrapper(self.wfield.func, fpool_l) 

400 for va, f in zip(fpool_l, F): 

401 vt = tuple(va) 

402 self[vt].f = f # set vertex object attribute v.f = f 

403 self.nfev += 1 

404 # Clean the pool 

405 self.fpool = set() 

406 

407 def pproc_fpool_nog(self): 

408 """ 

409 Process all field functions with no constraints supplied in parallel. 

410 """ 

411 self.wfield.func 

412 fpool_l = [] 

413 for v in self.fpool: 

414 fpool_l.append(v.x_a) 

415 F = self._mapwrapper(self.wfield.func, fpool_l) 

416 for va, f in zip(fpool_l, F): 

417 vt = tuple(va) 

418 self[vt].f = f # set vertex object attribute v.f = f 

419 self.nfev += 1 

420 # Clean the pool 

421 self.fpool = set() 

422 

423 def proc_minimisers(self): 

424 """Check for minimisers.""" 

425 for v in self: 

426 v.minimiser() 

427 v.maximiser() 

428 

429 

430class ConstraintWrapper: 

431 """Object to wrap constraints to pass to `multiprocessing.Pool`.""" 

432 def __init__(self, g_cons, g_cons_args): 

433 self.g_cons = g_cons 

434 self.g_cons_args = g_cons_args 

435 

436 def gcons(self, v_x_a): 

437 vfeasible = True 

438 for g, args in zip(self.g_cons, self.g_cons_args): 

439 # constraint may return more than 1 value. 

440 if np.any(g(v_x_a, *args) < 0.0): 

441 vfeasible = False 

442 break 

443 return vfeasible 

444 

445 

446class FieldWrapper: 

447 """Object to wrap field to pass to `multiprocessing.Pool`.""" 

448 def __init__(self, field, field_args): 

449 self.field = field 

450 self.field_args = field_args 

451 

452 def func(self, v_x_a): 

453 try: 

454 v_f = self.field(v_x_a, *self.field_args) 

455 except Exception: 

456 v_f = np.inf 

457 if np.isnan(v_f): 

458 v_f = np.inf 

459 

460 return v_f