193ca29b6SJames Wright // SPDX-FileCopyrightText: Copyright (c) 2017-2024, HONEE contributors. 293ca29b6SJames Wright // SPDX-License-Identifier: Apache-2.0 OR BSD-2-Clause 393ca29b6SJames Wright 493ca29b6SJames Wright /// @file 593ca29b6SJames Wright /// Custom file I/O functions for HONEE 693ca29b6SJames Wright 793ca29b6SJames Wright #include <honee-file.h> 893ca29b6SJames Wright 9c9ff4f08SJames Wright const PetscInt32 HONEE_FILE_TOKEN = 0xceedf00; // for backwards compatibility 10c9ff4f08SJames Wright const PetscInt32 HONEE_FILE_TOKEN_32 = 0xceedf32; 11c9ff4f08SJames Wright const PetscInt32 HONEE_FILE_TOKEN_64 = 0xceedf64; 1293ca29b6SJames Wright 13c9ff4f08SJames Wright // @brief Read in binary int based on it's data type 1493ca29b6SJames Wright static PetscErrorCode BinaryReadIntoInt(PetscViewer viewer, PetscInt *out, PetscDataType file_type) { 1593ca29b6SJames Wright PetscFunctionBeginUser; 1693ca29b6SJames Wright *out = -13; // appease the overzealous GCC compiler warning Gods 1793ca29b6SJames Wright if (file_type == PETSC_INT32) { 1893ca29b6SJames Wright PetscInt32 val; 1993ca29b6SJames Wright PetscCall(PetscViewerBinaryRead(viewer, &val, 1, NULL, PETSC_INT32)); 2093ca29b6SJames Wright *out = val; 2193ca29b6SJames Wright } else if (file_type == PETSC_INT64) { 2293ca29b6SJames Wright PetscInt64 val; 2393ca29b6SJames Wright PetscCall(PetscViewerBinaryRead(viewer, &val, 1, NULL, PETSC_INT64)); 2493ca29b6SJames Wright *out = val; 2593ca29b6SJames Wright } else { 2693ca29b6SJames Wright PetscCall(PetscViewerBinaryRead(viewer, out, 1, NULL, PETSC_INT)); 2793ca29b6SJames Wright } 2893ca29b6SJames Wright PetscFunctionReturn(PETSC_SUCCESS); 2993ca29b6SJames Wright } 3093ca29b6SJames Wright 31c9ff4f08SJames Wright /** 32c9ff4f08SJames Wright @brief Load vector from binary file, possibly with embedded solution time and step number 33c9ff4f08SJames Wright 34c9ff4f08SJames Wright Reads in Vec from binary file, possibly written by HONEE. 35c9ff4f08SJames Wright If written by HONEE, will also load the solution time and timestep, otherwise not. 36c9ff4f08SJames Wright Also handles case where file was written with different PetscInt size than is read with. 37c9ff4f08SJames Wright 38*360b37c9SJames Wright @param[in] viewer `PetscViewer` to read the vec from. Must be a binary viewer 39c9ff4f08SJames Wright @param[out] Q `Vec` to read the data into 40c9ff4f08SJames Wright @param[out] time Solution time the Vec was written at, or `NULL`. Set to 0 if legacy file format. 41c9ff4f08SJames Wright @param[out] step_number Timestep number the Vec was written at, or `NULL`. Set to 0 if legacy file format. 42c9ff4f08SJames Wright **/ 43*360b37c9SJames Wright PetscErrorCode HoneeLoadBinaryVec(PetscViewer viewer, Vec Q, PetscReal *time, PetscInt *step_number) { 4493ca29b6SJames Wright PetscInt file_step_number; 4593ca29b6SJames Wright PetscInt32 token; 4693ca29b6SJames Wright PetscReal file_time; 4793ca29b6SJames Wright PetscDataType file_type = PETSC_INT32; 48*360b37c9SJames Wright MPI_Comm comm = PetscObjectComm((PetscObject)viewer); 4993ca29b6SJames Wright 5093ca29b6SJames Wright PetscFunctionBeginUser; 5193ca29b6SJames Wright PetscCall(PetscViewerBinaryRead(viewer, &token, 1, NULL, PETSC_INT32)); 52c9ff4f08SJames Wright if (token == HONEE_FILE_TOKEN_32 || token == HONEE_FILE_TOKEN_64 || 53c9ff4f08SJames Wright token == HONEE_FILE_TOKEN) { // New style format; we're reading a file with step number and time in the header 54c9ff4f08SJames Wright if (token == HONEE_FILE_TOKEN_32) file_type = PETSC_INT32; 55c9ff4f08SJames Wright else if (token == HONEE_FILE_TOKEN_64) file_type = PETSC_INT64; 5693ca29b6SJames Wright PetscCall(BinaryReadIntoInt(viewer, &file_step_number, file_type)); 5793ca29b6SJames Wright PetscCall(PetscViewerBinaryRead(viewer, &file_time, 1, NULL, PETSC_REAL)); 5893ca29b6SJames Wright } else if (token == VEC_FILE_CLASSID) { // Legacy format of just the vector, encoded as [VEC_FILE_CLASSID, length, ] 5993ca29b6SJames Wright PetscInt length, N; 6093ca29b6SJames Wright PetscCall(BinaryReadIntoInt(viewer, &length, file_type)); 6193ca29b6SJames Wright PetscCall(VecGetSize(Q, &N)); 6293ca29b6SJames Wright PetscCheck(length == N, comm, PETSC_ERR_ARG_INCOMP, "File Vec has length %" PetscInt_FMT " but DM has global Vec size %" PetscInt_FMT, length, N); 6393ca29b6SJames Wright PetscCall(PetscViewerBinarySetSkipHeader(viewer, PETSC_TRUE)); 64c9ff4f08SJames Wright file_time = 0; 65c9ff4f08SJames Wright file_step_number = 0; 6693ca29b6SJames Wright } else SETERRQ(comm, PETSC_ERR_FILE_UNEXPECTED, "Not a fluids header token or a PETSc Vec in file"); 6793ca29b6SJames Wright 68c9ff4f08SJames Wright if (time) *time = file_time; 69c9ff4f08SJames Wright if (step_number) *step_number = file_step_number; 7093ca29b6SJames Wright PetscCall(VecLoad(Q, viewer)); 7193ca29b6SJames Wright PetscFunctionReturn(PETSC_SUCCESS); 7293ca29b6SJames Wright } 7393ca29b6SJames Wright 74c9ff4f08SJames Wright /** 75c9ff4f08SJames Wright @brief Open a PHASTA *.dat file, grabbing dimensions and file pointer 76c9ff4f08SJames Wright 77c9ff4f08SJames Wright This function opens the file specified by `path` using `PetscFOpen` and passes the file pointer in `fp`. 78c9ff4f08SJames Wright It is not closed in this function, thus `fp` must be closed sometime after this function has been called (using `PetscFClose` for example). 79c9ff4f08SJames Wright 80c9ff4f08SJames Wright Assumes that the first line of the file has the number of rows and columns as the only two entries, separated by a single space. 81c9ff4f08SJames Wright 82c9ff4f08SJames Wright @param[in] comm MPI_Comm for the program 83c9ff4f08SJames Wright @param[in] path Path to the file 84c9ff4f08SJames Wright @param[in] char_array_len Length of the character array that should contain each line 85c9ff4f08SJames Wright @param[out] dims Dimensions of the file, taken from the first line of the file 86c9ff4f08SJames Wright @param[out] fp File pointer to the opened file 87c9ff4f08SJames Wright **/ 8893ca29b6SJames Wright PetscErrorCode PhastaDatFileOpen(const MPI_Comm comm, const char path[PETSC_MAX_PATH_LEN], const PetscInt char_array_len, PetscInt dims[2], 8993ca29b6SJames Wright FILE **fp) { 9093ca29b6SJames Wright int ndims; 9193ca29b6SJames Wright char line[char_array_len]; 9293ca29b6SJames Wright char **array; 9393ca29b6SJames Wright 9493ca29b6SJames Wright PetscFunctionBeginUser; 9593ca29b6SJames Wright PetscCall(PetscFOpen(comm, path, "r", fp)); 9693ca29b6SJames Wright PetscCall(PetscSynchronizedFGets(comm, *fp, char_array_len, line)); 9793ca29b6SJames Wright PetscCall(PetscStrToArray(line, ' ', &ndims, &array)); 9893ca29b6SJames Wright PetscCheck(ndims == 2, comm, PETSC_ERR_FILE_UNEXPECTED, "Found %d dimensions instead of 2 on the first line of %s", ndims, path); 9993ca29b6SJames Wright 10093ca29b6SJames Wright for (PetscInt i = 0; i < ndims; i++) dims[i] = atoi(array[i]); 10193ca29b6SJames Wright PetscCall(PetscStrToArrayDestroy(ndims, array)); 10293ca29b6SJames Wright PetscFunctionReturn(PETSC_SUCCESS); 10393ca29b6SJames Wright } 10493ca29b6SJames Wright 105c9ff4f08SJames Wright /** 106c9ff4f08SJames Wright @brief Get the number of rows for the PHASTA file at path. 107c9ff4f08SJames Wright 108c9ff4f08SJames Wright Assumes that the first line of the file has the number of rows and columns as the only two entries, separated by a single space. 109c9ff4f08SJames Wright 110c9ff4f08SJames Wright @param[in] comm `MPI_Comm` for the program 111c9ff4f08SJames Wright @param[in] path Path to the file 112c9ff4f08SJames Wright @param[out] nrows Number of rows 113c9ff4f08SJames Wright **/ 11493ca29b6SJames Wright PetscErrorCode PhastaDatFileGetNRows(const MPI_Comm comm, const char path[PETSC_MAX_PATH_LEN], PetscInt *nrows) { 11593ca29b6SJames Wright const PetscInt char_array_len = 512; 11693ca29b6SJames Wright PetscInt dims[2]; 11793ca29b6SJames Wright FILE *fp; 11893ca29b6SJames Wright 11993ca29b6SJames Wright PetscFunctionBeginUser; 12093ca29b6SJames Wright PetscCall(PhastaDatFileOpen(comm, path, char_array_len, dims, &fp)); 12193ca29b6SJames Wright *nrows = dims[0]; 12293ca29b6SJames Wright PetscCall(PetscFClose(comm, fp)); 12393ca29b6SJames Wright PetscFunctionReturn(PETSC_SUCCESS); 12493ca29b6SJames Wright } 12593ca29b6SJames Wright 126c9ff4f08SJames Wright /** 127c9ff4f08SJames Wright @brief Read PetscReal values from PHASTA file 128c9ff4f08SJames Wright 129c9ff4f08SJames Wright @param[in] comm `MPI_Comm` for the program 130c9ff4f08SJames Wright @param[in] path Path to the file 131c9ff4f08SJames Wright @param[out] array Pointer to allocated array of correct size 132c9ff4f08SJames Wright **/ 13393ca29b6SJames Wright PetscErrorCode PhastaDatFileReadToArrayReal(MPI_Comm comm, const char path[PETSC_MAX_PATH_LEN], PetscReal array[]) { 13493ca29b6SJames Wright PetscInt dims[2]; 13593ca29b6SJames Wright FILE *fp; 13693ca29b6SJames Wright const PetscInt char_array_len = 512; 13793ca29b6SJames Wright char line[char_array_len]; 13893ca29b6SJames Wright 13993ca29b6SJames Wright PetscFunctionBeginUser; 14093ca29b6SJames Wright PetscCall(PhastaDatFileOpen(comm, path, char_array_len, dims, &fp)); 14193ca29b6SJames Wright 14293ca29b6SJames Wright for (PetscInt i = 0; i < dims[0]; i++) { 14393ca29b6SJames Wright int ndims; 14493ca29b6SJames Wright char **row_array; 14593ca29b6SJames Wright 14693ca29b6SJames Wright PetscCall(PetscSynchronizedFGets(comm, fp, char_array_len, line)); 14793ca29b6SJames Wright PetscCall(PetscStrToArray(line, ' ', &ndims, &row_array)); 14893ca29b6SJames Wright PetscCheck(ndims == dims[1], comm, PETSC_ERR_FILE_UNEXPECTED, 14993ca29b6SJames Wright "Line %" PetscInt_FMT " of %s does not contain enough columns (%d instead of %" PetscInt_FMT ")", i, path, ndims, dims[1]); 15093ca29b6SJames Wright 15193ca29b6SJames Wright for (PetscInt j = 0; j < dims[1]; j++) array[i * dims[1] + j] = (PetscReal)atof(row_array[j]); 15293ca29b6SJames Wright PetscCall(PetscStrToArrayDestroy(ndims, row_array)); 15393ca29b6SJames Wright } 15493ca29b6SJames Wright 15593ca29b6SJames Wright PetscCall(PetscFClose(comm, fp)); 15693ca29b6SJames Wright PetscFunctionReturn(PETSC_SUCCESS); 15793ca29b6SJames Wright } 158