Coverage for pySDC/helpers/fieldsIO.py: 72%

1#!/usr/bin/env python3

2# -*- coding: utf-8 -*-

3"""

4Generic utility class to write and read cartesian grid field solutions into binary files.

5It implements the base file handler class :class:`FieldsIO`, that is specialized into :

7- :class:`Scalar` : for 0D fields (scalar) with a given number of variables

8- :class:`Rectilinear` : for fields on N-dimensional rectilinear grids

10While each file handler need to be setup with specific parameters (grid, ...),

11each written file can be read using the same interface implemented in the

12base abstract class.

14Example

15-------

16>>> import numpy as np

17>>> from pySDC.helpers.fieldsIO import Rectilinear

18>>>

19>>> # Write some fields in files

20>>> x = np.linspace(0, 1, 128)

21>>> y = np.linspace(0, 1, 64)

22>>> fOut = Rectilinear(np.float64, "file.pysdc")

23>>> fOut.setHeader(nVar=2, coords=[x, y])

24>>> fOut.initialize()

25>>> times = [0, 1, 2]

26>>> xGrid, yGrid = np.meshgrid(x, y, indexing="ij")

27>>> u0 = np.array([-1, 1]).reshape((-1, 1, 1))*xGrid*yGrid

28>>> # u0 has shape [2, nX, nY]

29>>> for t in times:

30>>> fOut.addField(t, t*u0)

31>>>

32>>> # Read the file using the generic interface

33>>> from pySDC.helpers.fieldsIO import FieldsIO

34>>> fIn = FieldsIO.fromFile("file.pysdc")

35>>> times = fIn.times

36>>> assert len(times) == fIn.nFields

37>>> tEnd, uEnd = fIn.readField(-1)

38>>> assert tEnd == times[-1]

40Notes

41-----

42🚀 :class:`Rectilinear` is compatible with a MPI-based cartesian decomposition.

43See :class:`pySDC.helpers.fieldsIO.writeFields_MPI` for an illustrative example.

45Warning

46-------

47To use MPI collective writing, you need to call first the class methods :class:`Rectilinear.setupMPI` (cf their docstring).

48Also, `Rectilinear.setHeader` **must be given the global grids coordinates**, whether the code is run in parallel or not.

49"""

50import os

51import numpy as np

52from typing import Type, TypeVar

53import logging

55T = TypeVar("T")

57try:

58 try:

59 import dolfin as df # noqa: F841 (for some reason, dolfin always needs to be imported before mpi4py)

60 except ImportError:

61 pass

62 from mpi4py import MPI

63 from mpi4py.util.dtlib import from_numpy_dtype as MPI_DTYPE

64except ImportError:

66 class MPI:

67 COMM_WORLD = None

68 Intracomm = T

69 File = T

70 Datatype = T

72 def MPI_DTYPE():

73 pass

76# Supported data types

77DTYPES = {

78 0: np.float64, # double precision

79 1: np.complex128,

80}

81try:

82 DTYPES.update(

83 {

84 2: np.float128, # quadruple precision

85 3: np.complex256,

86 }

87 )

88except AttributeError:

89 logging.getLogger('FieldsIO').debug('Warning: Quadruple precision not available on this machine')

90try:

91 DTYPES.update(

92 {

93 4: np.float32, # single precision

94 5: np.complex64,

95 }

96 )

97except AttributeError:

98 logging.getLogger('FieldsIO').debug('Warning: Single precision not available on this machine')

100DTYPES_AVAIL = {val: key for key, val in DTYPES.items()}

101

102# Header dtype

103H_DTYPE = np.int8

104T_DTYPE = np.float64

105

106

107class FieldsIO:

108 """Abstract IO file handler"""

109

110 STRUCTS = {} # supported structures, modified dynamically

111 sID = None # structure ID of the FieldsIO class, modified dynamically

112

113 tSize = T_DTYPE().itemsize

114

115 ALLOW_OVERWRITE = False

116

117 def __init__(self, dtype, fileName):

118 """

119 Parameters

120 ----------

121 dtype : np.dtype

122 The data type of the fields values.

123 fileName : str

124 File.

125 """

126 assert dtype in DTYPES_AVAIL, f"{dtype=} not available. Supported on this machine: {list(DTYPES_AVAIL.keys())}"

127 self.dtype = dtype

128 self.fileName = fileName

129 self.initialized = False

130

131 # Initialized by the setHeader abstract method

132 self.header = None

133 self.nItems = None # number of values (dof) stored into one field

134

135 def __str__(self):

136 return f"FieldsIO[{self.__class__.__name__}|{self.dtype.__name__}|file:{self.fileName}]<{hex(id(self))}>"

137

138 def __repr__(self):

139 return self.__str__()

140

141 @classmethod

142 def fromFile(cls, fileName):

143 """

144 Read a file storing fields, and return the `FieldsIO` of the appropriate

145 field type (structure).

146

147 Parameters

148 ----------

149 fileName : str

150 Name of the binary file.

151

152 Returns

153 -------

154 fieldsIO : :class:`FieldsIO`

155 The specialized `FieldsIO` adapted to the file.

156 """

157 assert os.path.isfile(fileName), f"not a file ({fileName})"

158 with open(fileName, "rb") as f:

159 STRUCT, DTYPE = np.fromfile(f, dtype=H_DTYPE, count=2)

160 fieldsIO: FieldsIO = cls.STRUCTS[STRUCT](DTYPES[DTYPE], fileName)

161 fieldsIO.readHeader(f)

162 fieldsIO.initialized = True

163 return fieldsIO

164

165 @property

166 def hBase(self) -> np.ndarray:

167 """Base header into numpy array format"""

168 return np.array([self.sID, DTYPES_AVAIL[self.dtype]], dtype=H_DTYPE)

169

170 @classmethod

171 def register(cls, sID):

172 """

173 Decorator used to register a new class FieldsIO specialized class

174

175 Parameters

176 ----------

177 sID : int

178 Unique identifyer for the file, used in the binary file.

179 Since it's encoded on a 8-bytes signed integer,

180 it must be between -128 and 127

181

182 Example

183 -------

184 >>> # New specialized FieldsIO class

185 >>> @FieldsIO.register(sID=31)

186 >>> class HexaMesh2D(FieldsIO):

187 >>> pass # ... implementation

188 """

189

190 def wrapper(registered: Type[T]) -> Type[T]:

191 assert (

192 sID not in cls.STRUCTS

193 ), f"struct ID already taken by {cls.STRUCTS[sID]}, cannot use it for {registered}"

194 cls.STRUCTS[sID] = registered

195 registered.sID = sID

196 return registered

197

198 return wrapper

199

200 def initialize(self):

201 """Initialize the file handler : create the file with header, removing any existing file with the same name"""

202 assert self.header is not None, "header must be set before initializing FieldsIO"

203 assert not self.initialized, "FieldsIO already initialized"

204

205 if not self.ALLOW_OVERWRITE:

206 if os.path.isfile(self.fileName):

207 raise FileExistsError(

208 f"file {self.fileName!r} already exists, use FieldsIO.ALLOW_OVERWRITE = True to allow overwriting"

209 )

210

211 with open(self.fileName, "w+b") as f:

212 self.hBase.tofile(f)

213 for array in self.hInfos:

214 array.tofile(f)

215 self.initialized = True

216

217 def setHeader(self, **params):

218 """(Abstract) Set the header before creating a new file to store the fields"""

219 raise NotImplementedError()

220

221 @property

222 def hInfos(self) -> list[np.ndarray]:

223 """(Abstract) Array representing the grid structure to be written in the binary file."""

224 raise NotImplementedError()

225

226 def readHeader(self, f):

227 """

228 (Abstract) Read the header from the file storing the fields.

229

230 Parameters

231 ----------

232 f : `_io.TextIOWrapper`

233 File to read the header from.

234 """

235 raise NotImplementedError()

236

237 @property

238 def hSize(self):

239 """Size of the full header (in bytes)"""

240 return self.hBase.nbytes + sum(hInfo.nbytes for hInfo in self.hInfos)

241

242 @property

243 def itemSize(self):

244 """Size of one field value (in bytes)"""

245 return self.dtype().itemsize

246

247 @property

248 def fSize(self):

249 """Full size of a field (in bytes)"""

250 return self.nItems * self.itemSize

251

252 @property

253 def fileSize(self):

254 """Current size of the file (in bytes)"""

255 return os.path.getsize(self.fileName)

256

257 def addField(self, time, field):

258 """

259 Append one field solution at the end of the file with one given time.

260

261 Parameters

262 ----------

263 time : float-like

264 The associated time of the field solution.

265 field : np.ndarray

266 The field values.

267 """

268 assert self.initialized, "cannot add field to a non initialized FieldsIO"

269 field = np.asarray(field)

270 assert field.dtype == self.dtype, f"expected {self.dtype} dtype, got {field.dtype}"

271 assert field.size == self.nItems, f"expected {self.nItems} values, got {field.size}"

272 with open(self.fileName, "ab") as f:

273 np.array(time, dtype=T_DTYPE).tofile(f)

274 field.tofile(f)

275

276 @property

277 def nFields(self):

278 """Number of fields currently stored in the binary file"""

279 return int((self.fileSize - self.hSize) // (self.tSize + self.fSize))

280

281 def formatIndex(self, idx):

282 """Utility method to format a fields index to a positional integer (negative starts from last field index, like python lists)"""

283 nFields = self.nFields

284 if idx < 0:

285 idx = nFields + idx

286 assert idx < nFields, f"cannot read index {idx} from {nFields} fields"

287 assert idx >= 0, f"cannot read index {idx-nFields} from {nFields} fields"

288 return idx

289

290 @property

291 def times(self):

292 """Vector of all times stored in the binary file"""

293 times = []

294 with open(self.fileName, "rb") as f:

295 f.seek(self.hSize)

296 for i in range(self.nFields):

297 t = np.fromfile(f, dtype=T_DTYPE, count=1, offset=0 if i == 0 else self.fSize)[0]

298 times.append(float(t))

299 return times

300

301 def time(self, idx):

302 """Time stored at a given field index"""

303 idx = self.formatIndex(idx)

304 offset = self.hSize + idx * (self.tSize + self.fSize)

305 with open(self.fileName, "rb") as f:

306 t = np.fromfile(f, dtype=T_DTYPE, count=1, offset=offset)[0]

307 return float(t)

308

309 def readField(self, idx):

310 """

311 Read one field stored in the binary file, corresponding to the given

312 time index.

313

314 Parameters

315 ----------

316 idx : int

317 Positional index of the field.

318

319 Returns

320 -------

321 t : float

322 Stored time for this field.

323 field : np.ndarray

324 Read fields in a numpy array.

325 """

326 idx = self.formatIndex(idx)

327 offset = self.hSize + idx * (self.tSize + self.fSize)

328 with open(self.fileName, "rb") as f:

329 f.seek(offset)

330 t = float(np.fromfile(f, dtype=T_DTYPE, count=1)[0])

331 field = np.fromfile(f, dtype=self.dtype, count=self.nItems)

332 self.reshape(field)

333 return t, field

334

335 def reshape(self, field):

336 """Eventually reshape the field to correspond to the grid structure"""

337 pass

338

339

340@FieldsIO.register(sID=0)

341class Scalar(FieldsIO):

342 """FieldsIO handler storing a given number of scalar"""

343

344 # -------------------------------------------------------------------------

345 # Overridden methods

346 # -------------------------------------------------------------------------

347 def setHeader(self, nVar):

348 """

349 Set the descriptive grid structure to be stored in the file header.

350

351 Parameters

352 ----------

353 nVar : int

354 Number of scalar variable stored.

355 """

356 self.header = {"nVar": int(nVar)}

357 self.nItems = self.nVar

358

359 @property

360 def hInfos(self):

361 """Array representing the grid structure to be written in the binary file."""

362 return [np.array([self.nVar], dtype=np.int64)]

363

364 def readHeader(self, f):

365 """

366 Read the header from the binary file storing the fields.

367

368 Parameters

369 ----------

370 f : `_io.TextIOWrapper`

371 File to read the header from.

372 """

373 (nVar,) = np.fromfile(f, dtype=np.int64, count=1)

374 self.setHeader(nVar)

375

376 # -------------------------------------------------------------------------

377 # Class specifics

378 # -------------------------------------------------------------------------

379 @property

380 def nVar(self):

381 """Number of variables in a fields, as described in the header"""

382 return self.header["nVar"]

383

384

385@FieldsIO.register(sID=1)

386class Rectilinear(Scalar):

387 """FieldsIO handler storing a given number of scalar variables on a N-dimensional rectilinear grid"""

388

389 @staticmethod

390 def setupCoords(*coords):

391 """Utility function to setup grids in multiple dimensions, given the keyword arguments"""

392 coords = [np.asarray(coord, dtype=np.float64) for coord in coords]

393 for axis, coord in enumerate(coords):

394 assert coord.ndim == 1, f"coord for {axis=} must be one dimensional"

395 return coords

396

397 # -------------------------------------------------------------------------

398 # Overridden methods

399 # -------------------------------------------------------------------------

400 def setHeader(self, nVar, coords):

401 """

402 Set the descriptive grid structure to be stored in the file header.

403

404 Parameters

405 ----------

406 nVar : int

407 Number of 1D variables stored.

408 coords : np.1darray or list[np.1darray]

409 The grid coordinates in each dimensions.

410

411 Note

412 ----

413 When used in MPI decomposition mode, all coordinate **must** be the global grid.

414 """

415 if not isinstance(coords, (tuple, list)):

416 coords = [coords]

417 coords = self.setupCoords(*coords)

418 self.header = {"nVar": int(nVar), "coords": coords}

419 self.nItems = nVar * self.nDoF

420

421 @property

422 def hInfos(self):

423 """Array representing the grid structure to be written in the binary file."""

424 return [np.array([self.nVar, self.dim, *self.gridSizes], dtype=np.int32)] + [

425 np.array(coord, dtype=np.float64) for coord in self.header["coords"]

426 ]

427

428 def readHeader(self, f):

429 """

430 Read the header from the binary file storing the fields.

431

432 Parameters

433 ----------

434 f : `_io.TextIOWrapper`

435 File to read the header from.

436 """

437 nVar, dim = np.fromfile(f, dtype=np.int32, count=2)

438 gridSizes = np.fromfile(f, dtype=np.int32, count=dim)

439 coords = [np.fromfile(f, dtype=np.float64, count=n) for n in gridSizes]

440 self.setHeader(nVar, coords)

441

442 def reshape(self, fields: np.ndarray):

443 """Reshape the fields to a N-d array (inplace operation)"""

444 fields.shape = (self.nVar, *self.gridSizes)

445

446 # -------------------------------------------------------------------------

447 # Class specifics

448 # -------------------------------------------------------------------------

449 @property

450 def gridSizes(self):

451 """Number of points in y direction"""

452 return [coord.size for coord in self.header["coords"]]

453

454 @property

455 def dim(self):

456 """Number of grid dimensions"""

457 return len(self.gridSizes)

458

459 @property

460 def nDoF(self):

461 """Number of degrees of freedom for one variable"""

462 return np.prod(self.gridSizes)

463

464 def toVTR(self, baseName, varNames, idxFormat="{:06d}"):

465 """

466 Convert all 3D fields stored in binary format (FieldsIO) into a list

467 of VTR files, that can be read later with Paraview or equivalent to

468 make videos.

469

470 Parameters

471 ----------

472 baseName : str

473 Base name of the VTR file.

474 varNames : list[str]

475 Variable names of the fields.

476 idxFormat : str, optional

477 Formatting string for the index of the VTR file.

478 The default is "{:06d}".

479

480 Example

481 -------

482 >>> # Suppose the FieldsIO object is already written into outputs.pysdc

483 >>> import os

484 >>> from pySDC.utils.fieldsIO import Rectilinear

485 >>> os.makedirs("vtrFiles") # to store all VTR files into a subfolder

486 >>> Rectilinear.fromFile("outputs.pysdc").toVTR(

487 >>> baseName="vtrFiles/field", varNames=["u", "v", "w", "T", "p"])

488 """

489 assert self.dim == 3, "can only be used with 3D fields"

490 from pySDC.helpers.vtkIO import writeToVTR

491

492 template = f"{baseName}_{idxFormat}"

493 for i in range(self.nFields):

494 _, u = self.readField(i)

495 writeToVTR(template.format(i), u, self.header["coords"], varNames)

496

497 # -------------------------------------------------------------------------

498 # MPI-parallel implementation

499 # -------------------------------------------------------------------------

500 comm: MPI.Intracomm = None

501

502 @classmethod

503 def setupMPI(cls, comm: MPI.Intracomm, iLoc, nLoc):

504 """

505 Setup the MPI mode for the files IO, considering a decomposition

506 of the 1D grid into contiguous subintervals.

507

508 Parameters

509 ----------

510 comm : MPI.Intracomm

511 The space decomposition communicator.

512 iLoc : list[int]

513 Starting index of the local sub-domain in the global coordinates.

514 nLoc : list[int]

515 Number of points in the local sub-domain.

516 """

517 cls.comm = comm

518 cls.iLoc = iLoc

519 cls.nLoc = nLoc

520 cls.mpiFile: MPI.File = None

521 cls.mpiType: MPI.Datatype = None

522 cls.mpiFileType: MPI.Datatype = None

523

524 @property

525 def MPI_ON(self):

526 """Wether or not MPI is activated"""

527 if self.comm is None:

528 return False

529 return self.comm.Get_size() > 1

530

531 @property

532 def MPI_ROOT(self):

533 """Wether or not the process is MPI Root"""

534 if self.comm is None:

535 return True

536 return self.comm.Get_rank() == 0

537

538 def MPI_SETUP_FILETYPE(self):

539 """Setup subarray masks for each processes"""

540 self.mpiType = MPI_DTYPE(self.dtype)

541 self.mpiFileType = self.mpiType.Create_subarray(

542 [self.nVar, *self.gridSizes], # Global array sizes

543 [self.nVar, *self.nLoc], # Local array sizes

544 [0, *self.iLoc], # Global starting indices of local blocks

545 )

546 self.mpiFileType.Commit()

547

548 def MPI_FILE_OPEN(self, mode):

549 """Open the binary file in MPI mode"""

550 amode = {

551 "r": MPI.MODE_RDONLY,

552 "a": MPI.MODE_WRONLY | MPI.MODE_APPEND,

553 }[mode]

554 self.mpiFile = MPI.File.Open(self.comm, self.fileName, amode)

555 if self.mpiType is None:

556 self.MPI_SETUP_FILETYPE()

557

558 def MPI_WRITE(self, data):

559 """Write data (np.ndarray) in the binary file in MPI mode, at the current file cursor position."""

560 self.mpiFile.Write(data)

561

562 def MPI_WRITE_AT_ALL(self, offset, data: np.ndarray):

563 """

564 Write data in the binary file in MPI mode, with a given offset

565 **relative to the beginning of the file**.

566

567 Parameters

568 ----------

569 offset : int

570 Offset to write at, relative to the beginning of the file, in bytes.

571 data : np.ndarray

572 Data to be written in the binary file.

573 """

574 self.mpiFile.Set_view(disp=offset, etype=self.mpiType, filetype=self.mpiFileType)

575 self.mpiFile.Write_all(data)

576

577 def MPI_READ_AT_ALL(self, offset, data: np.ndarray):

578 """

579 Read data from the binary file in MPI mode, with a given offset

580 **relative to the beginning of the file**.

581

582 Parameters

583 ----------

584 offset : int

585 Offset to read at, relative to the beginning of the file, in bytes.

586 data : np.ndarray

587 Array on which to read the data from the binary file.

588 """

589 self.mpiFile.Set_view(disp=offset, etype=self.mpiType, filetype=self.mpiFileType)

590 self.mpiFile.Read_all(data)

591

592 def MPI_FILE_CLOSE(self):

593 """Close the binary file in MPI mode"""

594 self.mpiFile.Close()

595 self.mpiFile = None

596

597 def initialize(self):

598 """Initialize the binary file (write header) in MPI mode"""

599 if self.MPI_ROOT:

600 try:

601 super().initialize()

602 except AssertionError as e:

603 if self.MPI_ON:

604 print(f"{type(e)}: {e}")

605 self.comm.Abort()

606 else:

607 raise e

608

609 if self.MPI_ON:

610 self.comm.Barrier() # Important, should not be removed !

611 self.initialized = True

612

613 def addField(self, time, field):

614 """

615 Append one field solution at the end of the file with one given time,

616 possibly using MPI.

617

618 Parameters

619 ----------

620 time : float-like

621 The associated time of the field solution.

622 field : np.ndarray

623 The (local) field values.

624

625 Note

626 ----

627 If a MPI decomposition is used, field **must be** the local field values.

628 """

629 if not self.MPI_ON:

630 return super().addField(time, field)

631

632 assert self.initialized, "cannot add field to a non initialized FieldsIO"

633

634 field = np.asarray(field)

635 assert field.dtype == self.dtype, f"expected {self.dtype} dtype, got {field.dtype}"

636 assert field.shape == (

637 self.nVar,

638 *self.nLoc,

639 ), f"expected {(self.nVar, *self.nLoc)} shape, got {field.shape}"

640

641 offset = self.fileSize

642 self.MPI_FILE_OPEN(mode="a")

643

644 if self.MPI_ROOT:

645 self.MPI_WRITE(np.array(time, dtype=T_DTYPE))

646 offset += self.tSize

647 self.MPI_WRITE_AT_ALL(offset, field)

648 self.MPI_FILE_CLOSE()

649

650 def readField(self, idx):

651 """

652 Read one field stored in the binary file, corresponding to the given

653 time index, using MPI in the eventuality of space parallel decomposition.

654

655 Parameters

656 ----------

657 idx : int

658 Positional index of the field.

659

660 Returns

661 -------

662 t : float

663 Stored time for this field.

664 field : np.ndarray

665 Read (local) fields in a numpy array.

666

667 Note

668 ----

669 If a MPI decomposition is used, it reads and returns the local fields values only.

670 """

671 if not self.MPI_ON:

672 return super().readField(idx)

673

674 idx = self.formatIndex(idx)

675 offset = self.hSize + idx * (self.tSize + self.fSize)

676 with open(self.fileName, "rb") as f:

677 t = float(np.fromfile(f, dtype=T_DTYPE, count=1, offset=offset)[0])

678 offset += self.tSize

679

680 field = np.empty((self.nVar, *self.nLoc), dtype=self.dtype)

681

682 self.MPI_FILE_OPEN(mode="r")

683 self.MPI_READ_AT_ALL(offset, field)

684 self.MPI_FILE_CLOSE()

685

686 return t, field

687

688

689# -----------------------------------------------------------------------------------------------

690# Utility functions used for testing

691# -----------------------------------------------------------------------------------------------

692def initGrid(nVar, gridSizes):

693 dim = len(gridSizes)

694 coords = [np.linspace(0, 1, num=n, endpoint=False) for n in gridSizes]

695 s = [None] * dim

696 u0 = np.array(np.arange(nVar) + 1)[(slice(None), *s)]

697 for x in np.meshgrid(*coords, indexing="ij"):

698 u0 = u0 * x

699 return coords, u0

700

701

702def writeFields_MPI(fileName, dtypeIdx, algo, nSteps, nVar, gridSizes):

703 coords, u0 = initGrid(nVar, gridSizes)

704

705 from mpi4py import MPI

706 from pySDC.helpers.blocks import BlockDecomposition

707 from pySDC.helpers.fieldsIO import Rectilinear

708

709 comm = MPI.COMM_WORLD

710 MPI_SIZE = comm.Get_size()

711 MPI_RANK = comm.Get_rank()

712

713 blocks = BlockDecomposition(MPI_SIZE, gridSizes, algo, MPI_RANK)

714

715 iLoc, nLoc = blocks.localBounds

716 Rectilinear.setupMPI(comm, iLoc, nLoc)

717 s = [slice(i, i + n) for i, n in zip(iLoc, nLoc)]

718 u0 = u0[(slice(None), *s)]

719

720 f1 = Rectilinear(DTYPES[dtypeIdx], fileName)

721 f1.setHeader(nVar=nVar, coords=coords)

722

723 u0 = np.asarray(u0, dtype=f1.dtype)

724 f1.initialize()

725

726 times = np.arange(nSteps) / nSteps

727 for t in times:

728 ut = (u0 * t).astype(f1.dtype)

729 f1.addField(t, ut)

730

731 return u0

732

733

734def compareFields_MPI(fileName, u0, nSteps):

735 from pySDC.helpers.fieldsIO import FieldsIO

736

737 comm = MPI.COMM_WORLD

738 MPI_RANK = comm.Get_rank()

739 if MPI_RANK == 0:

740 print("Comparing fields with MPI")

741

742 f2 = FieldsIO.fromFile(fileName)

743

744 times = np.arange(nSteps) / nSteps

745 for idx, t in enumerate(times):

746 u1 = u0 * t

747 t2, u2 = f2.readField(idx)

748 assert t2 == t, f"fields[{idx}] in {f2} has incorrect time ({t2} instead of {t})"

749 assert u2.shape == u1.shape, f"{idx}'s fields in {f2} has incorrect shape"

750 assert np.allclose(u2, u1), f"{idx}'s fields in {f2} has incorrect values"