diff --git a/.github/workflows/mdf_verification.yml b/.github/workflows/mdf_verification.yml index a3c38126..a394a1e0 100644 --- a/.github/workflows/mdf_verification.yml +++ b/.github/workflows/mdf_verification.yml @@ -28,13 +28,14 @@ jobs: id: setup-python - name: Install Python dependencies - run: pip install --upgrade build pytest Bokeh numpy + run: + pip install --upgrade build pytest Bokeh numpy matplotlib - name: Download OpenFAST shell: bash run: | cd ${{github.workspace}}/ - git clone --recursive https://github.com/OpenFAST/openfast.git + git clone --recursive --single-branch -b dev https://github.com/OpenFAST/openfast.git - name: Compile MoorDynF shell: bash @@ -73,3 +74,11 @@ jobs: - name: Run the tests run: python ${{github.workspace}}/tests/.mdf_verification/verify.py ${{github.workspace}}/ + + - name: Upload the artifacts + uses: actions/upload-artifact@v4 + if: failure() + with: + name: MD-F comparison plots + path: | + ${{github.workspace}}/*.png \ No newline at end of file diff --git a/.gitignore b/.gitignore index 3be253f2..d0426d04 100644 --- a/.gitignore +++ b/.gitignore @@ -47,6 +47,13 @@ compile/SO/test* #program output files *.out +*.png +*.log +*.dat +tests/.mdf_verification/*.out +tests/.mdf_verification/*.log +tests/.mdf_verification/*.png +tests/.mdf_verification/*.dat # Python packaging # #################### diff --git a/docs/inputs.rst b/docs/inputs.rst index d977a1f4..318dc06c 100644 --- a/docs/inputs.rst +++ b/docs/inputs.rst @@ -595,6 +595,7 @@ The list of possible options is: If this is enabled initial conditions are calculated with scaled drag according to CdScaleIC. The new stationary solver in MoorDyn-C is more stable and more precise than the dynamic solver, but it can take longer to reach equilibrium. + - print_time (1): MoorDyn-C switch for printing the current timestep to the console. A note about time steps in MoorDyn-C: The internal time step is first taken from the dtM option. If no CFL factor is provided, then the user provided time step is used to calculate CFL and MoorDyn-C @@ -629,6 +630,7 @@ The following MoorDyn-C options are not supported by MoorDyn-F: - FricDamp: Same as CV in MoorDyn-F. - StatDynFricScale: Same as MC in MoorDyn-F. - ICgenDynamic: MoorDyn-F does not have a stationary solver for initial conditions + - print_time: Console printing in MoorDyn-F is controled by OpenFAST The following options from MoorDyn-F are not supported by MoorDyn-C: diff --git a/source/MoorDyn2.cpp b/source/MoorDyn2.cpp index a6a89a37..7478ea89 100644 --- a/source/MoorDyn2.cpp +++ b/source/MoorDyn2.cpp @@ -91,6 +91,7 @@ moordyn::MoorDyn::MoorDyn(const char* infilename, int log_level) , dtOut(0.0) , _t_integrator(NULL) , ICgenDynamic(false) + , print_time(true) , env(std::make_shared()) , GroundBody(NULL) , waves(nullptr) @@ -611,12 +612,15 @@ moordyn::MoorDyn::Step(const double* x, double& dt) { // should check if wave kinematics have been set up if expected! - const auto default_precision{std::cout.precision()}; - std::cout << std::fixed << setprecision(1); - LOGDBG << "t = " << t << "s \r"; - std::cout << std::defaultfloat << setprecision(default_precision); + + if (print_time) { + const auto default_precision{std::cout.precision()}; + std::cout << std::fixed << setprecision(1); + LOGDBG << "t = " << t << "s \r"; + std::cout << std::defaultfloat << setprecision(default_precision); - cout << "\rt = " << t << " " << flush; + cout << "\rt = " << t << " " << flush; + } if (dt <= 0) { // Nothing to do, just recover the forces if there are coupled DOFs @@ -2183,6 +2187,8 @@ moordyn::MoorDyn::readOptionsLine(vector& in_txt, int i) this->seafloor->setup(env, filepath); } else if (name == "ICgenDynamic") ICgenDynamic = bool(atof(entries[0].c_str())); + else if (name == "print_time") + print_time = bool(atof(entries[0].c_str())); else LOGWRN << "Warning: Unrecognized option '" << name << "'" << endl; } diff --git a/source/MoorDyn2.hpp b/source/MoorDyn2.hpp index 599a8234..b2ba79d7 100644 --- a/source/MoorDyn2.hpp +++ b/source/MoorDyn2.hpp @@ -550,6 +550,8 @@ class MoorDyn final : public io::IO real ICthresh; // use dynamic (true) or stationary (false) inital condition solver bool ICgenDynamic; + /// print the timesteps. Disable for MD_F reg tests + bool print_time; // temporary wave kinematics flag used to store input value while keeping // env.WaveKin=0 for IC gen moordyn::waves::waves_settings WaveKinTemp; diff --git a/tests/.mdf_verification/verify.py b/tests/.mdf_verification/verify.py index 14c39ef9..4d548368 100644 --- a/tests/.mdf_verification/verify.py +++ b/tests/.mdf_verification/verify.py @@ -85,17 +85,19 @@ def to_num(s): if line.startswith('"'): # The value is a string - end = line[1:].find('"') + 1 - value = line[1:end] + end = line[1:].find('"') + 2 + value = line[1:end-1] else: - end = line.find(' ') + indicies = [line.find(' '), line.find(' '), line.find('\t')] + if -1 in indicies: + indicies.remove(-1) + end = min(indicies) value = to_num(line[:end]) - - line = line[end:] + line = line[end:].strip() while line.find(" ") != -1: line = line.replace(" ", " ") - key = line.split(" ")[1] + key = line.split(" ")[0] data[key] = value return data @@ -106,8 +108,11 @@ def read_driver(test): test_root = os.path.join(args.root, 'openfast/reg_tests/r-test/modules/moordyn/', test) - with open(os.path.join(test_root, "md_driver.inp"), "r") as f: - lines = f.readlines() + try: + with open(os.path.join(test_root, "md_driver.inp"), "r") as f: + lines = f.readlines() + except FileNotFoundError: + return -1,-1 def get_section(lines, name): start, end = None, None @@ -153,7 +158,7 @@ def create_input_file(env, md): end = i break opts = mdopts2dict(opts[:end]) - for optin, optout in (('Gravity', 'gravity'), + for optin, optout in (('Gravity', 'gravity'), # todo: make case insensitive ('rhoW', 'rho'), ('WtrDpth', 'WtrDpth')): if optin not in env.keys(): @@ -168,6 +173,8 @@ def create_input_file(env, md): f"{optout}={opts[optout]} in MoorDyn config file") lines.insert(start + end, f"{env[optin]} {optout} option set by the driver\n") + lines.insert(start+end, "0 print_time - console printing of timesteps disabled for GH actions" ) + for line in lines: fout.write(line) return fname @@ -254,34 +261,79 @@ def to_num(s): return float(s) data = [] + heads = [] with open(fpath, "r") as fin: - lines = fin.readlines()[skiplines:] - for line in lines: + lines = fin.readlines()[skiplines-2:] + line = lines[0] + line = line.strip().replace("\t", " ") + while line.find(" ") != -1: + line = line.replace(" ", " ") + heads = [field for field in line.split()] + + for line in lines[2:]: line = line.strip().replace("\t", " ") while line.find(" ") != -1: line = line.replace(" ", " ") data.append([to_num(field) for field in line.split()]) - return np.transpose(data) + return np.transpose(data), heads -def plot(ref, data, fpath): +def plot(ref, ref_heads, data, data_heads, passed, rtol_mag, atol_mag, fpath): if plt is None: return - colors = list(mcolors.XKCD_COLORS.values()) + fig,axes = plt.subplots(ref.shape[0]-1, 2, sharex = True, figsize=(12.8,4*(ref.shape[0]-1))) + if ref.shape[0] <= 2: + ax = [axes] + else: + ax = axes for i in range(1, ref.shape[0]): - plt.plot(ref[0, :], ref[i, :], linestyle='dashed', - color=colors[i - 1]) - plt.plot(data[0, :], data[i, :], linestyle='solid', - color=colors[i - 1], label=f'channel {i}') - plt.legend(loc='best') - plt.savefig(fpath) + # Plot the channels + ax[i-1][0].plot(data[0, :], data[i, :], linestyle='solid', + color='b', label=f"MD-C: {(data_heads[i])}") + ax[i-1][0].plot(ref[0, :], ref[i, :], linestyle='dashed', + color='r', label=f"MD-F: {(ref_heads[i])}") + ax[i-1][0].legend(loc='best') + + # Plot the difference. First calculate threshold (from OpenFAST testing erorPlotting.py _plotError) + NUMEPS = 1e-12 + ATOL_MIN = 1e-6 + baseline_offset = ref[i, :] - np.min(ref[i, :]) + b_order_of_magnitude = np.floor( np.log10( baseline_offset + NUMEPS ) ) + rtol = 10**(-1 * rtol_mag) + atol = 10**(max(b_order_of_magnitude) - atol_mag) + atol = max(atol, ATOL_MIN) + passfail_line = atol + rtol * abs(ref[i, :]) + ax[i-1][1].plot(data[0, :], passfail_line, color='g', + label=f"Threshold {(data_heads[i])}") + ax[i-1][1].plot(data[0, :], abs(ref[i, :] - data[i, :]), color='gray', + label=f"Error {(data_heads[i])}") + ax[i-1][1].legend(loc='best') + + # titles + if i == 1: + if passed: + ax[i-1][0].set_title("PASSED") + else: + ax[i-1][0].set_title("FAILED") + + ax[i-1][1].set_title(f"abs diff btwn MD-C and MD-F") + + + + fig.tight_layout() + fig.savefig(fpath, dpi=400) + plt.close() # Run the tests... summary = {} for test in tests: + print("\n ------------------") print(f"Test {test}...") env, md = read_driver(test) + if env == md == -1: + print(f"WARNING: Error in loading {test}. Does not match md_driver.inp convention") + continue fname = create_input_file(env, md) system = moordyn.Create(fname) # Get the NDoFs and check if the motions are right @@ -305,17 +357,18 @@ def plot(ref, data, fpath): moordyn.Step(system, rorg, u, t, dt) moordyn.Close(system) # Read the ouputs and compare - ref = read_outs(md["OutRootName"], skiplines=8) - new = read_outs(os.path.splitext(fname)[0] + ".out", skiplines=2) - # Drop the eventual points at the tail that ight come from precision errors + ref, ref_heads = read_outs(md["OutRootName"], skiplines=8) + new, new_heads = read_outs(os.path.splitext(fname)[0] + ".out", skiplines=2) + # Drop the eventual points at the tail that might come from precision errors # on the time n_samples = min(ref.shape[1], new.shape[1]) ref = ref[:, :n_samples] new = new[:, :n_samples] - plot(ref, new, test + ".png") passing = np.all( pass_fail.passing_channels(ref, new, args.rtol, args.atol)) + plot(ref, ref_heads, new, new_heads, passing, args.rtol, args.atol, test + ".png") + summary[test] = passing print("")