From e39d42140c783df5f9e6e7a46e9ec73e58f437e2 Mon Sep 17 00:00:00 2001 From: Vikram Date: Tue, 28 Jun 2022 19:04:30 +0530 Subject: [PATCH] update code after rebase with master --- Project.toml | 1 + docs/pages.jl | 78 ++++++----- src/GEKPLS.jl | 336 +++++++++++++++++++++++----------------------- src/Surrogates.jl | 16 +-- test/GEKPLS.jl | 148 ++++++++++---------- test/runtests.jl | 22 +-- 6 files changed, 303 insertions(+), 298 deletions(-) diff --git a/Project.toml b/Project.toml index d87569e97..9893178ad 100644 --- a/Project.toml +++ b/Project.toml @@ -8,6 +8,7 @@ Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f" ExtendableSparse = "95c220a8-a1cf-11e9-0c77-dbfce5f500b3" GLM = "38e38edf-8417-5370-95a0-9cbb8c7f171a" IterativeSolvers = "42fd0dbc-a981-5370-80f2-aaf504508153" +JuliaFormatter = "98e50ef6-434e-11e9-1051-2b60c6c9e899" LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" QuasiMonteCarlo = "8a4e6c94-4038-4cdc-81c3-7e6ffdb2a71b" SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" diff --git a/docs/pages.jl b/docs/pages.jl index 5101a1e11..b53c2d048 100644 --- a/docs/pages.jl +++ b/docs/pages.jl @@ -1,40 +1,38 @@ -pages = [ - "index.md" - "Tutorials" => [ - "Basics" => "tutorials.md", - "Radials" => "radials.md", - "Kriging" => "kriging.md", - "Gaussian Process" =>"abstractgps.md", - "Lobachevsky" => "lobachevsky.md", - "Linear" => "LinearSurrogate.md", - "InverseDistance" => "InverseDistance.md", - "RandomForest" => "randomforest.md", - "SecondOrderPolynomial" => "secondorderpoly.md", - "NeuralSurrogate" => "neural.md", - "Wendland" => "wendland.md", - "Polynomial Chaos" => "polychaos.md", - "Variable Fidelity" => "variablefidelity.md", - "Gradient Enhanced Kriging" => "gek.md", - "GEKPLS" => "gekpls.md" - ] - "User guide" => [ - "Samples" => "samples.md", - "Surrogates" => "surrogate.md", - "Optimization" => "optimizations.md" - ] - "Benchmarks" => [ - "Sphere function" => "sphere_function.md", - "Lp norm" => "lp.md", - "Rosenbrock" => "rosenbrock.md", - "Tensor product" => "tensor_prod.md", - "Cantilever beam" => "cantilever.md", - "Water Flow function" => "water_flow.md", - "Welded beam function" => "welded_beam.md", - "Branin function" => "BraninFunction.md", - "Ackley function" => "ackley.md", - "Gramacy & Lee Function" => "gramacylee.md", - "Salustowicz Benchmark" => "Salustowicz.md", - "Multi objective optimization" => "multi_objective_opt.md" - ] - "Contributing" => "contributing.md" - ] +pages = ["index.md" + "Tutorials" => [ + "Basics" => "tutorials.md", + "Radials" => "radials.md", + "Kriging" => "kriging.md", + "Gaussian Process" => "abstractgps.md", + "Lobachevsky" => "lobachevsky.md", + "Linear" => "LinearSurrogate.md", + "InverseDistance" => "InverseDistance.md", + "RandomForest" => "randomforest.md", + "SecondOrderPolynomial" => "secondorderpoly.md", + "NeuralSurrogate" => "neural.md", + "Wendland" => "wendland.md", + "Polynomial Chaos" => "polychaos.md", + "Variable Fidelity" => "variablefidelity.md", + "Gradient Enhanced Kriging" => "gek.md", + "GEKPLS" => "gekpls.md", + ] + "User guide" => [ + "Samples" => "samples.md", + "Surrogates" => "surrogate.md", + "Optimization" => "optimizations.md", + ] + "Benchmarks" => [ + "Sphere function" => "sphere_function.md", + "Lp norm" => "lp.md", + "Rosenbrock" => "rosenbrock.md", + "Tensor product" => "tensor_prod.md", + "Cantilever beam" => "cantilever.md", + "Water Flow function" => "water_flow.md", + "Welded beam function" => "welded_beam.md", + "Branin function" => "BraninFunction.md", + "Ackley function" => "ackley.md", + "Gramacy & Lee Function" => "gramacylee.md", + "Salustowicz Benchmark" => "Salustowicz.md", + "Multi objective optimization" => "multi_objective_opt.md", + ] + "Contributing" => "contributing.md"] diff --git a/src/GEKPLS.jl b/src/GEKPLS.jl index c228dea9b..255cbc588 100644 --- a/src/GEKPLS.jl +++ b/src/GEKPLS.jl @@ -1,18 +1,18 @@ using LinearAlgebra using Statistics -mutable struct GEKPLS{T<:AbstractFloat} <: AbstractSurrogate +mutable struct GEKPLS{T <: AbstractFloat} <: AbstractSurrogate x::Matrix{T} #1 y::Matrix{T} #2 grads::Matrix{T} #3 xl::Matrix{T} #xlimits #4 - delta:: T #5 + delta::T #5 extra_points::Int #6 num_components::Int #7 beta::Vector{T} #8 gamma::Matrix{T} #9 theta::Vector{T} #10 - reduced_likelihood_function_value:: T #11 + reduced_likelihood_function_value::T #11 X_offset::Matrix{T} #12 X_scale::Matrix{T} #13 X_after_std::Matrix{T} #14 - X after standardization @@ -22,11 +22,11 @@ mutable struct GEKPLS{T<:AbstractFloat} <: AbstractSurrogate end function bounds_error(x, xl) - num_x_rows = size(x,1) + num_x_rows = size(x, 1) num_dim = size(xl, 1) for i in 1:num_x_rows for j in 1:num_dim - if (x[i, j] < xl[j,1] || x[i,j] > xl[j,2]) + if (x[i, j] < xl[j, 1] || x[i, j] > xl[j, 2]) return true end end @@ -36,7 +36,7 @@ end #constructor for GEKPLS Struct function GEKPLS(X, y, grads, n_comp, delta_x, xlimits, extra_points, θ) - + #ensure that X values are within the upper and lower bounds if bounds_error(X, xlimits) println("X values outside bounds") @@ -44,34 +44,38 @@ function GEKPLS(X, y, grads, n_comp, delta_x, xlimits, extra_points, θ) end theta = [θ for i in 1:n_comp] - pls_mean, X_after_PLS, y_after_PLS = _ge_compute_pls(X, y, n_comp, grads, delta_x, xlimits, extra_points) - X_after_std, y_after_std, X_offset, y_mean, X_scale, y_std = standardization(X_after_PLS, y_after_PLS) + pls_mean, X_after_PLS, y_after_PLS = _ge_compute_pls(X, y, n_comp, grads, delta_x, + xlimits, extra_points) + X_after_std, y_after_std, X_offset, y_mean, X_scale, y_std = standardization(X_after_PLS, + y_after_PLS) D, ij = cross_distances(X_after_std) - pls_mean_reshaped = reshape(pls_mean, (size(X,2),n_comp)) + pls_mean_reshaped = reshape(pls_mean, (size(X, 2), n_comp)) d = componentwise_distance_PLS(D, "squar_exp", n_comp, pls_mean_reshaped) nt, nd = size(X_after_PLS) - beta, gamma, reduced_likelihood_function_value = _reduced_likelihood_function(theta, "squar_exp", d, nt, ij, y_after_std) - return GEKPLS(X, y, grads, xlimits, delta_x, extra_points, n_comp, beta, gamma, theta, reduced_likelihood_function_value, - X_offset, X_scale, X_after_std, pls_mean_reshaped, y_mean, y_std) + beta, gamma, reduced_likelihood_function_value = _reduced_likelihood_function(theta, + "squar_exp", + d, nt, ij, + y_after_std) + return GEKPLS(X, y, grads, xlimits, delta_x, extra_points, n_comp, beta, gamma, theta, + reduced_likelihood_function_value, + X_offset, X_scale, X_after_std, pls_mean_reshaped, y_mean, y_std) println("struct created") end - # predictor function (g::GEKPLS)(X_test) n_eval, n_features_x = size(X_test) X_cont = (X_test .- g.X_offset) ./ g.X_scale dx = differences(X_cont, g.X_after_std) pred_d = componentwise_distance_PLS(dx, "squar_exp", g.num_components, g.pls_mean) - nt = size(g.X_after_std,1) - r = transpose(reshape(squar_exp(g.theta, pred_d),(nt,n_eval))) - f = ones(n_eval,1) + nt = size(g.X_after_std, 1) + r = transpose(reshape(squar_exp(g.theta, pred_d), (nt, n_eval))) + f = ones(n_eval, 1) y_ = (f * g.beta) + (r * g.gamma) y = g.y_mean .+ g.y_std * y_ return y end - function add_point!(g::GEKPLS, new_x, new_y, new_grads) if new_x in g.x println("Adding a sample that already exists. Cannot build GEKPLS") @@ -80,20 +84,27 @@ function add_point!(g::GEKPLS, new_x, new_y, new_grads) if bounds_error(new_x, g.xl) println("x values outside bounds") - return + return end - g.x = vcat(g.x, new_x) g.y = vcat(g.y, new_y) g.grads = vcat(g.grads, new_grads) - pls_mean, X_after_PLS, y_after_PLS = _ge_compute_pls(g.x, g.y, g.num_components, g.grads, g.delta, g.xl, g.extra_points) - g.X_after_std, y_after_std, g.X_offset, g.y_mean, g.X_scale, g.y_std = standardization(X_after_PLS, y_after_PLS) + pls_mean, X_after_PLS, y_after_PLS = _ge_compute_pls(g.x, g.y, g.num_components, + g.grads, g.delta, g.xl, + g.extra_points) + g.X_after_std, y_after_std, g.X_offset, g.y_mean, g.X_scale, g.y_std = standardization(X_after_PLS, + y_after_PLS) D, ij = cross_distances(g.X_after_std) g.pls_mean = reshape(pls_mean, (size(g.x, 2), g.num_components)) d = componentwise_distance_PLS(D, "squar_exp", g.num_components, g.pls_mean) nt, nd = size(X_after_PLS) - g.beta, g.gamma, g.reduced_likelihood_function_value = _reduced_likelihood_function(g.theta, "squar_exp", d, nt, ij, y_after_std) + g.beta, g.gamma, g.reduced_likelihood_function_value = _reduced_likelihood_function(g.theta, + "squar_exp", + d, + nt, + ij, + y_after_std) end function _ge_compute_pls(X, y, n_comp, grads, delta_x, xlimits, extra_points) @@ -119,45 +130,43 @@ function _ge_compute_pls(X, y, n_comp, grads, delta_x, xlimits, extra_points) # and https://github.com/SMTorg/smt/blob/f124c01ffa78c04b80221dded278a20123dac742/smt/surrogate_models/gekpls.py#L48 nt, dim = size(X) - XX = zeros(0,dim) - yy = zeros(0,size(y)[2]) + XX = zeros(0, dim) + yy = zeros(0, size(y)[2]) coeff_pls = zeros((dim, n_comp, nt)) - + for i in 1:nt - if dim >= 3 - bb_vals = circshift(boxbehnken(dim, 1),1) + if dim >= 3 + bb_vals = circshift(boxbehnken(dim, 1), 1) else - bb_vals = [ - 0.0 0.0; #center - 1.0 0.0; #right - 0.0 1.0; #up - -1.0 0.0; #left - 0.0 -1.0; #down - 1.0 1.0; #right up - -1.0 1.0; #left up - -1.0 -1.0; #left down - 1.0 -1.0; #right down - ] + bb_vals = [0.0 0.0; #center + 1.0 0.0; #right + 0.0 1.0; #up + -1.0 0.0; #left + 0.0 -1.0; #down + 1.0 1.0; #right up + -1.0 1.0; #left up + -1.0 -1.0; #left down + 1.0 -1.0] end - _X = zeros((size(bb_vals)[1], dim)) - _y = zeros((size(bb_vals)[1], 1)) + _X = zeros((size(bb_vals)[1], dim)) + _y = zeros((size(bb_vals)[1], 1)) bb_vals = bb_vals .* (delta_x * (xlimits[:, 2] - xlimits[:, 1]))' #smt calls this sign. I've called it bb_vals - _X = X[i, :]' .+ bb_vals - bb_vals = bb_vals .* grads[i,:]' - _y = y[i, :] .+ sum(bb_vals, dims=2) - + _X = X[i, :]' .+ bb_vals + bb_vals = bb_vals .* grads[i, :]' + _y = y[i, :] .+ sum(bb_vals, dims = 2) + #_pls.fit(_X, _y) # relic from sklearn versiom; retained for future reference. #coeff_pls[:, :, i] = _pls.x_rotations_ #relic from sklearn versiom; retained for future reference. coeff_pls[:, :, i] = _modified_pls(_X, _y, n_comp) #_modified_pls returns the equivalent of SKLearn's _pls.x_rotations_ if extra_points != 0 - start_index = max(1, length(coeff_pls[:,1,i])-extra_points+1) - max_coeff = sortperm(broadcast(abs,coeff_pls[:,1,i]))[start_index:end] + start_index = max(1, length(coeff_pls[:, 1, i]) - extra_points + 1) + max_coeff = sortperm(broadcast(abs, coeff_pls[:, 1, i]))[start_index:end] for ii in max_coeff XX = [XX; transpose(X[i, :])] - XX[end, ii] += delta_x * (xlimits[ii,2]-xlimits[ii,1]) + XX[end, ii] += delta_x * (xlimits[ii, 2] - xlimits[ii, 1]) yy = [yy; y[i]] - yy[end] += grads[i,ii] * delta_x * (xlimits[ii,2]-xlimits[ii,1]) + yy[end] += grads[i, ii] * delta_x * (xlimits[ii, 2] - xlimits[ii, 1]) end end end @@ -166,89 +175,89 @@ function _ge_compute_pls(X, y, n_comp, grads, delta_x, xlimits, extra_points) y = [y; yy] end - pls_mean = mean(broadcast(abs,coeff_pls),dims=3) + pls_mean = mean(broadcast(abs, coeff_pls), dims = 3) return pls_mean, X, y end ######start of bbdesign###### +# +# Adapted from 'ExperimentalDesign.jl: Design of Experiments in Julia' +# https://github.com/phrb/ExperimentalDesign.jl - # - # Adapted from 'ExperimentalDesign.jl: Design of Experiments in Julia' - # https://github.com/phrb/ExperimentalDesign.jl +# MIT License - # MIT License +# ExperimentalDesign.jl: Design of Experiments in Julia +# Copyright (C) 2019 Pedro Bruel - # ExperimentalDesign.jl: Design of Experiments in Julia - # Copyright (C) 2019 Pedro Bruel +# Permission is hereby granted, free of charge, to any person obtaining a copy of +# this software and associated documentation files (the "Software"), to deal in +# the Software without restriction, including without limitation the rights to +# use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of +# the Software, and to permit persons to whom the Software is furnished to do so, +# subject to the following conditions: - # Permission is hereby granted, free of charge, to any person obtaining a copy of - # this software and associated documentation files (the "Software"), to deal in - # the Software without restriction, including without limitation the rights to - # use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of - # the Software, and to permit persons to whom the Software is furnished to do so, - # subject to the following conditions: +# The above copyright notice and this permission notice (including the next +# paragraph) shall be included in all copies or substantial portions of the +# Software. - # The above copyright notice and this permission notice (including the next - # paragraph) shall be included in all copies or substantial portions of the - # Software. +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR +# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS +# FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR +# COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER +# IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN +# CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. +# - # THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR - # IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS - # FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR - # COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER - # IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN - # CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. - # +function boxbehnken(matrix_size::Int) + boxbehnken(matrix_size, matrix_size) +end - function boxbehnken(matrix_size::Int) - boxbehnken(matrix_size, matrix_size) - end - - function boxbehnken(matrix_size::Int, center::Int) - @assert matrix_size>=3 - - A_fact = explicit_fullfactorial(Tuple([-1,1] for i = 1:2)) - - rows = floor(Int, (0.5*matrix_size*(matrix_size-1))*size(A_fact)[1]) - - A = zeros(rows, matrix_size) - - l = 0 - for i in 1:matrix_size-1 - for j in i+1:matrix_size - l = l +1 - A[max(0, (l - 1)*size(A_fact)[1])+1:l*size(A_fact)[1], i] = A_fact[:, 1] - A[max(0, (l - 1)*size(A_fact)[1])+1:l*size(A_fact)[1], j] = A_fact[:, 2] - end - end - - if center == matrix_size - if matrix_size <= 16 - points = [0, 0, 3, 3, 6, 6, 6, 8, 9, 10, 12, 12, 13, 14, 15, 16] - center = points[matrix_size] - end +function boxbehnken(matrix_size::Int, center::Int) + @assert matrix_size >= 3 + + A_fact = explicit_fullfactorial(Tuple([-1, 1] for i in 1:2)) + + rows = floor(Int, (0.5 * matrix_size * (matrix_size - 1)) * size(A_fact)[1]) + + A = zeros(rows, matrix_size) + + l = 0 + for i in 1:(matrix_size - 1) + for j in (i + 1):matrix_size + l = l + 1 + A[(max(0, (l - 1) * size(A_fact)[1]) + 1):(l * size(A_fact)[1]), i] = A_fact[:, + 1] + A[(max(0, (l - 1) * size(A_fact)[1]) + 1):(l * size(A_fact)[1]), j] = A_fact[:, + 2] end - - A = transpose(hcat(transpose(A), transpose(zeros(center, matrix_size)))) - end - - function explicit_fullfactorial(factors::Tuple) - explicit_fullfactorial(fullfactorial(factors)) - end - - function explicit_fullfactorial(iterator::Base.Iterators.ProductIterator) - hcat(vcat.(collect(iterator)...)...) end - - function fullfactorial(factors::Tuple) - Base.Iterators.product(factors...) + + if center == matrix_size + if matrix_size <= 16 + points = [0, 0, 3, 3, 6, 6, 6, 8, 9, 10, 12, 12, 13, 14, 15, 16] + center = points[matrix_size] + end end - + + A = transpose(hcat(transpose(A), transpose(zeros(center, matrix_size)))) +end + +function explicit_fullfactorial(factors::Tuple) + explicit_fullfactorial(fullfactorial(factors)) +end + +function explicit_fullfactorial(iterator::Base.Iterators.ProductIterator) + hcat(vcat.(collect(iterator)...)...) +end + +function fullfactorial(factors::Tuple) + Base.Iterators.product(factors...) +end ######end of bb design###### -function standardization(X,y) +function standardization(X, y) """ We substract the mean from each variable. Then, we divide the values of each variable by its standard deviation. @@ -280,14 +289,13 @@ function standardization(X,y) #Equivalent of https://github.com/SMTorg/smt/blob/4a4df255b9259965439120091007f9852f41523e/smt/utils/kriging_utils.py#L21 X_offset = mean(X, dims = 1) X_scale = std(X, dims = 1) - X_scale = map(x -> (x==0.0) ? x=1 : x=x, X_scale) #to prevent division by 0 below + X_scale = map(x -> (x == 0.0) ? x = 1 : x = x, X_scale) #to prevent division by 0 below y_mean = mean(y) y_std = std(y) - y_std = map(y -> (y==0) ? y=1 : y=y, y_std) #to prevent division by 0 below - X = (X.-X_offset) ./ X_scale + y_std = map(y -> (y == 0) ? y = 1 : y = y, y_std) #to prevent division by 0 below + X = (X .- X_offset) ./ X_scale y = (y .- y_mean) ./ y_std return X, y, X_offset, y_mean, X_scale, y_std - end function cross_distances(X) @@ -311,18 +319,17 @@ function cross_distances(X) """ # equivalent of https://github.com/SMTorg/smt/blob/4a4df255b9259965439120091007f9852f41523e/smt/utils/kriging_utils.py#L86 n_samples, n_features = size(X) - n_nonzero_cross_dist = ( n_samples * (n_samples - 1) ) ÷ 2 + n_nonzero_cross_dist = (n_samples * (n_samples - 1)) ÷ 2 ij = zeros((n_nonzero_cross_dist, 2)) D = zeros((n_nonzero_cross_dist, n_features)) ll_1 = 0 - - for k in 1:n_samples - 1 + + for k in 1:(n_samples - 1) ll_0 = ll_1 + 1 ll_1 = ll_0 + n_samples - k - 1 ij[ll_0:ll_1, 1] .= k - ij[ll_0:ll_1, 2] = k+1:1:n_samples - D[ll_0:ll_1, :] = -(X[(k + 1) : n_samples,:] .- X[k,:]') - + ij[ll_0:ll_1, 2] = (k + 1):1:n_samples + D[ll_0:ll_1, :] = -(X[(k + 1):n_samples, :] .- X[k, :]') end return D, Int.(ij) end @@ -367,7 +374,7 @@ function componentwise_distance_PLS(D, corr, n_comp, coeff_pls) D_corr = zeros((size(D)[1], n_comp)) if corr == "squar_exp" - D_corr = D.^2 * coeff_pls.^2 + D_corr = D .^ 2 * coeff_pls .^ 2 else #abs_exp D_corr = abs.(D) * abs.(coeff_pls) end @@ -387,27 +394,27 @@ function squar_exp(theta, d) Returns: -------- r: array containing the values of the autocorrelation model - + """ n_components = size(d)[2] - theta = reshape(theta, (1,n_components)) - return exp.(-sum(theta .* d, dims=2)) + theta = reshape(theta, (1, n_components)) + return exp.(-sum(theta .* d, dims = 2)) end function differences(X, Y) #equivalent of https://github.com/SMTorg/smt/blob/4a4df255b9259965439120091007f9852f41523e/smt/utils/kriging_utils.py#L392 #code credit: Elias Carvalho - https://stackoverflow.com/questions/72392010/row-wise-operations-between-matrices-in-julia - Rx = repeat(X, inner=(size(Y, 1), 1)) + Rx = repeat(X, inner = (size(Y, 1), 1)) Ry = repeat(Y, size(X, 1)) return Rx - Ry end -function _reduced_likelihood_function(theta, kernel_type, d, nt, ij, y_norma, noise=0.0) +function _reduced_likelihood_function(theta, kernel_type, d, nt, ij, y_norma, noise = 0.0) """ This function is a loose translation of SMT code from https://github.com/SMTorg/smt/blob/4a4df255b9259965439120091007f9852f41523e/smt/surrogate_models/krg_based.py#L247 It determines the BLUP parameters and evaluates the reduced likelihood function for the given theta. - + Parameters ---------- theta: array containing the parameters at which the Gaussian Process model parameters should be determined. @@ -430,68 +437,65 @@ function _reduced_likelihood_function(theta, kernel_type, d, nt, ij, y_norma, no #equivalent of https://github.com/SMTorg/smt/blob/4a4df255b9259965439120091007f9852f41523e/smt/surrogate_models/krg_based.py#L247 reduced_likelihood_function_value = -Inf nugget = 1000000.0 * eps() #a jitter for numerical stability; reducing the multiple from 1000000.0 results in positive definite error for Cholesky decomposition; - if kernel_type =="squar_exp" #todo - add other kernel type abs_exp etc. - r = squar_exp(theta,d) + if kernel_type == "squar_exp" #todo - add other kernel type abs_exp etc. + r = squar_exp(theta, d) end - R = (I + zeros(nt,nt)) .* (1.0 + nugget + noise) + R = (I + zeros(nt, nt)) .* (1.0 + nugget + noise) - for k in 1:size(ij)[1] - R[ij[k,1],ij[k,2]]=r[k] - R[ij[k,2],ij[k,1]]=r[k] + for k in 1:size(ij)[1] + R[ij[k, 1], ij[k, 2]] = r[k] + R[ij[k, 2], ij[k, 1]] = r[k] end C = cholesky(R).L #todo - #values diverge at this point from SMT code; verify impact - F = ones(nt,1) #todo - examine if this should be a parameter for this function - Ft = C\F + F = ones(nt, 1) #todo - examine if this should be a parameter for this function + Ft = C \ F Q, G = qr(Ft) Q = Array(Q) - Yt = C\y_norma + Yt = C \ y_norma #todo - in smt, they check if the matrix is ill-conditioned using SVD. Verify and include if necessary beta = G \ [(transpose(Q) ⋅ Yt)] rho = Yt .- (Ft .* beta) gamma = transpose(C) \ rho - sigma2 = sum((rho).^2, dims=1)/nt - detR = prod(diag(C).^(2.0/nt)) - reduced_likelihood_function_value = - nt * log10(sum(sigma2)) - nt * log10(detR) + sigma2 = sum((rho) .^ 2, dims = 1) / nt + detR = prod(diag(C) .^ (2.0 / nt)) + reduced_likelihood_function_value = -nt * log10(sum(sigma2)) - nt * log10(detR) return beta, gamma, reduced_likelihood_function_value end - ### MODIFIED PLS BELOW ### # The code below is a simplified version of # SKLearn's PLS # https://github.com/scikit-learn/scikit-learn/blob/80598905e/sklearn/cross_decomposition/_pls.py - -function _center_scale(X,Y) +function _center_scale(X, Y) x_mean = mean(X, dims = 1) X .-= x_mean y_mean = mean(Y, dims = 1) Y .-= y_mean - x_std = std(X, dims = 1) + x_std = std(X, dims = 1) x_std[x_std .== 0] .= 1.0 X ./= x_std - y_std = std(Y, dims = 1) - y_std[y_std .==0] .= 1.0 + y_std = std(Y, dims = 1) + y_std[y_std .== 0] .= 1.0 Y ./= y_std - return X,Y + return X, Y end -function _svd_flip_1d(u,v) +function _svd_flip_1d(u, v) # equivalent of https://github.com/scikit-learn/scikit-learn/blob/80598905e517759b4696c74ecc35c6e2eb508cff/sklearn/cross_decomposition/_pls.py#L149 biggest_abs_val_idx = findmax(abs.(vec(u)))[2] sign_ = sign(u[biggest_abs_val_idx]) u .*= sign_ v .*= sign_ - end -function _get_first_singular_vectors_power_method(X,Y) +function _get_first_singular_vectors_power_method(X, Y) my_eps = eps() y_score = vec(Y) x_weights = transpose(X)y_score / dot(y_score, y_score) - x_weights ./= (sqrt(dot(x_weights,x_weights)) + my_eps) + x_weights ./= (sqrt(dot(x_weights, x_weights)) + my_eps) x_score = X * x_weights y_weights = transpose(Y)x_score / dot(x_score, x_score) y_score = Y * y_weights / (dot(y_weights, y_weights) + my_eps) @@ -502,28 +506,28 @@ function _get_first_singular_vectors_power_method(X,Y) return x_weights, y_weights end -function _modified_pls(X,Y, n_components) - x_weights_ = zeros(size(X,2),n_components) - _x_scores = zeros(size(X,1), n_components) - x_loadings_ = zeros(size(X,2),n_components) - Xk, Yk = _center_scale(X,Y) +function _modified_pls(X, Y, n_components) + x_weights_ = zeros(size(X, 2), n_components) + _x_scores = zeros(size(X, 1), n_components) + x_loadings_ = zeros(size(X, 2), n_components) + Xk, Yk = _center_scale(X, Y) for k in 1:n_components - x_weights, y_weights = _get_first_singular_vectors_power_method(Xk,Yk) + x_weights, y_weights = _get_first_singular_vectors_power_method(Xk, Yk) if x_weights == false break - end + end - _svd_flip_1d(x_weights, y_weights) + _svd_flip_1d(x_weights, y_weights) x_scores = Xk * x_weights x_loadings = transpose(x_scores)Xk / dot(x_scores, x_scores) Xk = Xk - (x_scores * x_loadings) - y_loadings = transpose(x_scores)*Yk / dot(x_scores, x_scores) - Yk = Yk - x_scores*y_loadings - x_weights_[:,k] = x_weights - _x_scores[:,k] = x_scores - x_loadings_[:,k] = vec(x_loadings) + y_loadings = transpose(x_scores) * Yk / dot(x_scores, x_scores) + Yk = Yk - x_scores * y_loadings + x_weights_[:, k] = x_weights + _x_scores[:, k] = x_scores + x_loadings_[:, k] = vec(x_loadings) end x_rotations_ = x_weights_ * pinv(transpose(x_loadings_)x_weights_) diff --git a/src/Surrogates.jl b/src/Surrogates.jl index 16335f97a..eaf50ed44 100644 --- a/src/Surrogates.jl +++ b/src/Surrogates.jl @@ -31,11 +31,11 @@ function RadialBasisStructure(; radial_function, scale_factor, sparse) end #Kriging structure: -function KrigingStructure(;p,theta) +function KrigingStructure(; p, theta) return (name = "Kriging", p = p, theta = theta) end -function GEKStructure(;p,theta) +function GEKStructure(; p, theta) return (name = "GEK", p = p, theta = theta) end @@ -45,12 +45,12 @@ function LinearStructure() end #InverseDistance structure -function InverseDistanceStructure(;p) +function InverseDistanceStructure(; p) return (name = "InverseDistanceSurrogate", p = p) end #Lobachevsky structure -function LobachevskyStructure(;alpha,n,sparse) +function LobachevskyStructure(; alpha, n, sparse) return (name = "LobachevskySurrogate", alpha = alpha, n = n, sparse = sparse) end @@ -61,7 +61,7 @@ function NeuralStructure(; model, loss, opt, n_echos) end #Random forest structure -function RandomForestStructure(;num_round) +function RandomForestStructure(; num_round) return (name = "RandomForestSurrogate", num_round = num_round) end @@ -81,7 +81,7 @@ function PolyChaosStructure(; op) end export current_surrogates -export GEKPLS +export GEKPLS export RadialBasisStructure, KrigingStructure, LinearStructure, InverseDistanceStructure export LobachevskyStructure, NeuralStructure, RandomForestStructure, SecondOrderPolynomialStructure @@ -89,7 +89,7 @@ export WendlandStructure export AbstractSurrogate, SamplingAlgorithm export Kriging, RadialBasis, add_point!, current_estimate, std_error_at_point # radial basis functions -export linearRadial,cubicRadial,multiquadricRadial,thinplateRadial +export linearRadial, cubicRadial, multiquadricRadial, thinplateRadial # samplers export sample, GridSample, UniformSample, SobolSample, LatinHypercubeSample, @@ -97,7 +97,7 @@ export sample, GridSample, UniformSample, SobolSample, LatinHypercubeSample, export RandomSample, KroneckerSample, GoldenSample, SectionSample # Optimization algorithms -export SRBF,LCBS,EI,DYCORS,SOP,EGO,RTEA,SMB,surrogate_optimize +export SRBF, LCBS, EI, DYCORS, SOP, EGO, RTEA, SMB, surrogate_optimize export LobachevskySurrogate, lobachevsky_integral, lobachevsky_integrate_dimension export LinearSurrogate export SVMSurrogate diff --git a/test/GEKPLS.jl b/test/GEKPLS.jl index 28c070989..26f33f39b 100644 --- a/test/GEKPLS.jl +++ b/test/GEKPLS.jl @@ -1,7 +1,6 @@ using Surrogates using Zygote - function vector_of_tuples_to_matrix(v) #convert training data generated by surrogate sampling into a matrix suitable for GEKPLS num_rows = length(v) @@ -9,7 +8,7 @@ function vector_of_tuples_to_matrix(v) K = zeros(num_rows, num_cols) for row in 1:num_rows for col in 1:num_cols - K[row, col]=v[row][col] + K[row, col] = v[row][col] end end return K @@ -38,55 +37,56 @@ function water_flow(x) H_l = x[6] L = x[7] K_w = x[8] - log_val = log(r/r_w) - return (2*pi*T_u*(H_u - H_l))/ ( log_val*(1 + (2*L*T_u/(log_val*r_w^2*K_w)) + T_u/T_l)) + log_val = log(r / r_w) + return (2 * pi * T_u * (H_u - H_l)) / + (log_val * (1 + (2 * L * T_u / (log_val * r_w^2 * K_w)) + T_u / T_l)) end n = 1000 d = 8 -lb = [0.05,100,63070,990,63.1,700,1120,9855] -ub = [0.15,50000,115600,1110,116,820,1680,12045] -x = sample(n,lb,ub,SobolSample()) +lb = [0.05, 100, 63070, 990, 63.1, 700, 1120, 9855] +ub = [0.15, 50000, 115600, 1110, 116, 820, 1680, 12045] +x = sample(n, lb, ub, SobolSample()) X = vector_of_tuples_to_matrix(x) grads = vector_of_tuples_to_matrix2(gradient.(water_flow, x)) -y = reshape(water_flow.(x),(size(x,1),1)) +y = reshape(water_flow.(x), (size(x, 1), 1)) xlimits = hcat(lb, ub) -n_test = 100 -x_test = sample(n_test,lb,ub,GoldenSample()) -X_test = vector_of_tuples_to_matrix(x_test) +n_test = 100 +x_test = sample(n_test, lb, ub, GoldenSample()) +X_test = vector_of_tuples_to_matrix(x_test) y_true = water_flow.(x_test) -@testset "Test 1: Water Flow Function Test (dimensions = 8; n_comp = 2; extra_points = 2)" begin +@testset "Test 1: Water Flow Function Test (dimensions = 8; n_comp = 2; extra_points = 2)" begin n_comp = 2 delta_x = 0.0001 extra_points = 2 initial_theta = 0.01 - g = GEKPLS(X, y, grads, n_comp, delta_x, xlimits, extra_points, initial_theta) + g = GEKPLS(X, y, grads, n_comp, delta_x, xlimits, extra_points, initial_theta) y_pred = g(X_test) - rmse = sqrt(sum(((y_pred - y_true).^2)/n_test)) - @test isapprox(rmse, 0.03, atol=0.02) #rmse: 0.039 + rmse = sqrt(sum(((y_pred - y_true) .^ 2) / n_test)) + @test isapprox(rmse, 0.03, atol = 0.02) #rmse: 0.039 end -@testset "Test 2: Water Flow Function Test (dimensions = 8; n_comp = 3; extra_points = 2)" begin +@testset "Test 2: Water Flow Function Test (dimensions = 8; n_comp = 3; extra_points = 2)" begin n_comp = 3 delta_x = 0.0001 extra_points = 2 initial_theta = 0.01 g = GEKPLS(X, y, grads, n_comp, delta_x, xlimits, extra_points, initial_theta) #change hard-coded 2 param to variable y_pred = g(X_test) - rmse = sqrt(sum(((y_pred - y_true).^2)/n_test)) - @test isapprox(rmse, 0.02, atol=0.01) #rmse: 0.027 + rmse = sqrt(sum(((y_pred - y_true) .^ 2) / n_test)) + @test isapprox(rmse, 0.02, atol = 0.01) #rmse: 0.027 end -@testset "Test 3: Water Flow Function Test (dimensions = 8; n_comp = 3; extra_points = 3)" begin +@testset "Test 3: Water Flow Function Test (dimensions = 8; n_comp = 3; extra_points = 3)" begin n_comp = 3 delta_x = 0.0001 extra_points = 3 initial_theta = 0.01 - g = GEKPLS(X, y, grads, n_comp, delta_x, xlimits, extra_points, initial_theta) + g = GEKPLS(X, y, grads, n_comp, delta_x, xlimits, extra_points, initial_theta) y_pred = g(X_test) - rmse = sqrt(sum(((y_pred - y_true).^2)/n_test)) - @test isapprox(rmse, 0.02, atol=0.01) #rmse: 0.027 + rmse = sqrt(sum(((y_pred - y_true) .^ 2) / n_test)) + @test isapprox(rmse, 0.02, atol = 0.01) #rmse: 0.027 end ## welded beam tests @@ -94,76 +94,77 @@ function welded_beam(x) h = x[1] l = x[2] t = x[3] - a = 6000/(sqrt(2)*h*l) - b = (6000*(14+0.5*l)*sqrt(0.25*(l^2+(h+t)^2)))/(2*(0.707*h*l*(l^2/12 + 0.25*(h+t)^2))) - return (sqrt(a^2+b^2 + l*a*b))/(sqrt(0.25*(l^2+(h+t)^2))) + a = 6000 / (sqrt(2) * h * l) + b = (6000 * (14 + 0.5 * l) * sqrt(0.25 * (l^2 + (h + t)^2))) / + (2 * (0.707 * h * l * (l^2 / 12 + 0.25 * (h + t)^2))) + return (sqrt(a^2 + b^2 + l * a * b)) / (sqrt(0.25 * (l^2 + (h + t)^2))) end n = 1000 d = 3 -lb = [0.125,5.0,5.0] -ub = [1.,10.,10.] -x = sample(n,lb,ub,SobolSample()) +lb = [0.125, 5.0, 5.0] +ub = [1.0, 10.0, 10.0] +x = sample(n, lb, ub, SobolSample()) X = vector_of_tuples_to_matrix(x) grads = vector_of_tuples_to_matrix2(gradient.(welded_beam, x)) -y = reshape(welded_beam.(x),(size(x,1),1)) +y = reshape(welded_beam.(x), (size(x, 1), 1)) xlimits = hcat(lb, ub) -n_test = 100 -x_test = sample(n_test,lb,ub,GoldenSample()) -X_test = vector_of_tuples_to_matrix(x_test) +n_test = 100 +x_test = sample(n_test, lb, ub, GoldenSample()) +X_test = vector_of_tuples_to_matrix(x_test) y_true = welded_beam.(x_test) -@testset "Test 4: Welded Beam Function Test (dimensions = 3; n_comp = 3; extra_points = 2)" begin +@testset "Test 4: Welded Beam Function Test (dimensions = 3; n_comp = 3; extra_points = 2)" begin n_comp = 3 delta_x = 0.0001 extra_points = 2 initial_theta = 0.01 - g = GEKPLS(X, y, grads, n_comp, delta_x, xlimits, extra_points, initial_theta) + g = GEKPLS(X, y, grads, n_comp, delta_x, xlimits, extra_points, initial_theta) y_pred = g(X_test) - rmse = sqrt(sum(((y_pred - y_true).^2)/n_test)) - @test isapprox(rmse, 39.0, atol=0.5) #rmse: 38.988 + rmse = sqrt(sum(((y_pred - y_true) .^ 2) / n_test)) + @test isapprox(rmse, 39.0, atol = 0.5) #rmse: 38.988 end -@testset "Test 5: Welded Beam Function Test (dimensions = 3; n_comp = 2; extra_points = 2)" begin +@testset "Test 5: Welded Beam Function Test (dimensions = 3; n_comp = 2; extra_points = 2)" begin n_comp = 2 delta_x = 0.0001 extra_points = 2 initial_theta = 0.01 - g = GEKPLS(X, y, grads, n_comp, delta_x, xlimits, extra_points, initial_theta) + g = GEKPLS(X, y, grads, n_comp, delta_x, xlimits, extra_points, initial_theta) y_pred = g(X_test) - rmse = sqrt(sum(((y_pred - y_true).^2)/n_test)) - @test isapprox(rmse, 39.5, atol=0.5) #rmse: 39.481 + rmse = sqrt(sum(((y_pred - y_true) .^ 2) / n_test)) + @test isapprox(rmse, 39.5, atol = 0.5) #rmse: 39.481 end ## increasing extra points increases accuracy -@testset "Test 6: Welded Beam Function Test (dimensions = 3; n_comp = 2; extra_points = 4)" begin +@testset "Test 6: Welded Beam Function Test (dimensions = 3; n_comp = 2; extra_points = 4)" begin n_comp = 2 delta_x = 0.0001 extra_points = 4 initial_theta = 0.01 - g = GEKPLS(X, y, grads, n_comp, delta_x, xlimits, extra_points, initial_theta) + g = GEKPLS(X, y, grads, n_comp, delta_x, xlimits, extra_points, initial_theta) y_pred = g(X_test) - rmse = sqrt(sum(((y_pred - y_true).^2)/n_test)) - @test isapprox(rmse, 37.5, atol=0.5) #rmse: 37.87 + rmse = sqrt(sum(((y_pred - y_true) .^ 2) / n_test)) + @test isapprox(rmse, 37.5, atol = 0.5) #rmse: 37.87 end ## sphere function tests function sphere_function(x) - return sum(x.^2) + return sum(x .^ 2) end ## 3D n = 100 lb = [-5.0, -5.0, -5.0] -ub = [5.0, 5.0 ,5.0] -x = sample(n,lb,ub,SobolSample()) +ub = [5.0, 5.0, 5.0] +x = sample(n, lb, ub, SobolSample()) X = vector_of_tuples_to_matrix(x) grads = vector_of_tuples_to_matrix2(gradient.(sphere_function, x)) -y = reshape(sphere_function.(x),(size(x,1),1)) +y = reshape(sphere_function.(x), (size(x, 1), 1)) xlimits = hcat(lb, ub) -n_test = 100 -x_test = sample(n_test,lb,ub,GoldenSample()) -X_test = vector_of_tuples_to_matrix(x_test) +n_test = 100 +x_test = sample(n_test, lb, ub, GoldenSample()) +X_test = vector_of_tuples_to_matrix(x_test) y_true = sphere_function.(x_test) @testset "Test 7: Sphere Function Test (dimensions = 3; n_comp = 2; extra_points = 2)" begin @@ -171,25 +172,25 @@ y_true = sphere_function.(x_test) delta_x = 0.0001 extra_points = 2 initial_theta = 0.01 - g = GEKPLS(X, y, grads, n_comp, delta_x, xlimits, extra_points, initial_theta) + g = GEKPLS(X, y, grads, n_comp, delta_x, xlimits, extra_points, initial_theta) y_pred = g(X_test) - rmse = sqrt(sum(((y_pred - y_true).^2)/n_test)) - @test isapprox(rmse, 0.001, atol=0.05) #rmse: 0.00083 + rmse = sqrt(sum(((y_pred - y_true) .^ 2) / n_test)) + @test isapprox(rmse, 0.001, atol = 0.05) #rmse: 0.00083 end ## 2D -n = 50 +n = 50 d = 2 lb = [-10.0, -10.0] ub = [10.0, 10.0] -x = sample(n,lb,ub,SobolSample()) +x = sample(n, lb, ub, SobolSample()) X = vector_of_tuples_to_matrix(x) grads = vector_of_tuples_to_matrix2(gradient.(sphere_function, x)) -y = reshape(sphere_function.(x),(size(x,1),1)) +y = reshape(sphere_function.(x), (size(x, 1), 1)) xlimits = hcat(lb, ub) -n_test = 10 -x_test = sample(n_test,lb,ub,GoldenSample()) -X_test = vector_of_tuples_to_matrix(x_test) +n_test = 10 +x_test = sample(n_test, lb, ub, GoldenSample()) +X_test = vector_of_tuples_to_matrix(x_test) y_true = sphere_function.(x_test) @testset "Test 8: Sphere Function Test (dimensions = 2; n_comp = 2; extra_points = 2" begin @@ -197,16 +198,16 @@ y_true = sphere_function.(x_test) delta_x = 0.0001 extra_points = 2 initial_theta = 0.01 - g = GEKPLS(X, y, grads, n_comp, delta_x, xlimits, extra_points, initial_theta) + g = GEKPLS(X, y, grads, n_comp, delta_x, xlimits, extra_points, initial_theta) y_pred = g(X_test) - rmse = sqrt(sum(((y_pred - y_true).^2)/n_test)) - @test isapprox(rmse, 0.1, atol=0.5) #rmse: 0.0022 + rmse = sqrt(sum(((y_pred - y_true) .^ 2) / n_test)) + @test isapprox(rmse, 0.1, atol = 0.5) #rmse: 0.0022 end @testset "Test 9: Add Point Test (dimensions = 3; n_comp = 2; extra_points = 2)" begin #first we create a surrogate model with just 3 input points initial_x_vec = [(1.0, 2.0, 3.0), (4.0, 5.0, 6.0), (7.0, 8.0, 9.0)] - initial_y = reshape(sphere_function.(initial_x_vec), (size(initial_x_vec,1),1)) + initial_y = reshape(sphere_function.(initial_x_vec), (size(initial_x_vec, 1), 1)) initial_X = vector_of_tuples_to_matrix(initial_x_vec) initial_grads = vector_of_tuples_to_matrix2(gradient.(sphere_function, initial_x_vec)) lb = [-5.0, -5.0, -5.0] @@ -216,22 +217,23 @@ end delta_x = 0.0001 extra_points = 2 initial_theta = 0.01 - g = GEKPLS(initial_X, initial_y, initial_grads, n_comp, delta_x, xlimits, extra_points, initial_theta) - n_test = 100 - x_test = sample(n_test,lb,ub,GoldenSample()) - X_test = vector_of_tuples_to_matrix(x_test) - y_true = sphere_function.(x_test) + g = GEKPLS(initial_X, initial_y, initial_grads, n_comp, delta_x, xlimits, extra_points, + initial_theta) + n_test = 100 + x_test = sample(n_test, lb, ub, GoldenSample()) + X_test = vector_of_tuples_to_matrix(x_test) + y_true = sphere_function.(x_test) y_pred1 = g(X_test) - rmse1 = sqrt(sum(((y_pred1 - y_true).^2)/n_test)) #rmse1 = 31.91 + rmse1 = sqrt(sum(((y_pred1 - y_true) .^ 2) / n_test)) #rmse1 = 31.91 #then we update the model with more points to see if performance improves n = 100 - x = sample(n,lb,ub,SobolSample()) + x = sample(n, lb, ub, SobolSample()) X = vector_of_tuples_to_matrix(x) grads = vector_of_tuples_to_matrix2(gradient.(sphere_function, x)) - y = reshape(sphere_function.(x),(size(x,1),1)) + y = reshape(sphere_function.(x), (size(x, 1), 1)) add_point!(g, X, y, grads) y_pred2 = g(X_test) - rmse2 = sqrt(sum(((y_pred2 - y_true).^2)/n_test)) #rmse2 = 0.0015 + rmse2 = sqrt(sum(((y_pred2 - y_true) .^ 2) / n_test)) #rmse2 = 0.0015 @test (rmse2 < rmse1) end diff --git a/test/runtests.jl b/test/runtests.jl index aa703cca4..e435b486f 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -4,8 +4,8 @@ using SafeTestsets using Pkg function dev_subpkg(subpkg) - subpkg_path = joinpath(dirname(@__DIR__), "lib", subpkg) - Pkg.develop(PackageSpec(path=subpkg_path)) + subpkg_path = joinpath(dirname(@__DIR__), "lib", subpkg) + Pkg.develop(PackageSpec(path = subpkg_path)) end for pkg in ["SurrogatesAbstractGPs", "SurrogatesFlux", "SurrogatesPolyChaos", "SurrogatesRandomForest", "SurrogatesSVM"] @@ -15,15 +15,15 @@ for pkg in ["SurrogatesAbstractGPs", "SurrogatesFlux", "SurrogatesPolyChaos", end end -@time @safetestset "GEKPLS.jl" begin include("GEKPLS.jl") end -@time @safetestset "Radials.jl" begin include("Radials.jl") end -@time @safetestset "Kriging.jl" begin include("Kriging.jl") end -@time @safetestset "Sampling" begin include("sampling.jl") end -@time @safetestset "Optimization" begin include("optimization.jl") end -@time @safetestset "LinearSurrogate" begin include("linearSurrogate.jl") end -@time @safetestset "Lobachevsky" begin include("lobachevsky.jl") end -@time @safetestset "InverseDistanceSurrogate" begin include("inverseDistanceSurrogate.jl") end -@time @safetestset "SecondOrderPolynomialSurrogate" begin include("secondOrderPolynomialSurrogate.jl") end +@time @safetestset "GEKPLS.jl" begin include("GEKPLS.jl") end +@time @safetestset "Radials.jl" begin include("Radials.jl") end +@time @safetestset "Kriging.jl" begin include("Kriging.jl") end +@time @safetestset "Sampling" begin include("sampling.jl") end +@time @safetestset "Optimization" begin include("optimization.jl") end +@time @safetestset "LinearSurrogate" begin include("linearSurrogate.jl") end +@time @safetestset "Lobachevsky" begin include("lobachevsky.jl") end +@time @safetestset "InverseDistanceSurrogate" begin include("inverseDistanceSurrogate.jl") end +@time @safetestset "SecondOrderPolynomialSurrogate" begin include("secondOrderPolynomialSurrogate.jl") end # @time @safetestset "AD_Compatibility" begin include("AD_compatibility.jl") end @time @safetestset "Wendland" begin include("Wendland.jl") end @time @safetestset "VariableFidelity" begin include("VariableFidelity.jl") end