*.config
*.creator
+*.creator.user
*.files
*.includes
+++ /dev/null
-/*
- * This file is part of the GROMACS molecular simulation package.
- *
- * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
- * Copyright (c) 2001-2004, The GROMACS development team.
- * Copyright (c) 2013,2014,2015,2016, by the GROMACS development team, led by
- * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
- * and including many others, as listed in the AUTHORS file in the
- * top-level source directory and at http://www.gromacs.org.
- *
- * GROMACS is free software; you can redistribute it and/or
- * modify it under the terms of the GNU Lesser General Public License
- * as published by the Free Software Foundation; either version 2.1
- * of the License, or (at your option) any later version.
- *
- * GROMACS is distributed in the hope that it will be useful,
- * but WITHOUT ANY WARRANTY; without even the implied warranty of
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
- * Lesser General Public License for more details.
- *
- * You should have received a copy of the GNU Lesser General Public
- * License along with GROMACS; if not, see
- * http://www.gnu.org/licenses, or write to the Free Software Foundation,
- * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
- *
- * If you want to redistribute modifications to GROMACS, please
- * consider that scientific software is very special. Version
- * control is crucial - bugs must be traceable. We will be happy to
- * consider code for inclusion in the official distribution, but
- * derived work must not be called official GROMACS. Details are found
- * in the README & COPYING files - if they are missing, get the
- * official version at http://www.gromacs.org.
- *
- * To help us fund GROMACS development, we humbly ask that you cite
- * the research papers on the package. Check out http://www.gromacs.org.
- */
-#ifndef GMX_MATH_VECTYPES_H
-#define GMX_MATH_VECTYPES_H
-
-#include <math.h> // since this file used from C source files TODO: replace with cmath
-
-#include "gromacs/utility/real.h"
-
-#define XX 0 /* Defines for indexing in */
-#define YY 1 /* vectors */
-#define ZZ 2
-#define DIM 3 /* Dimension of vectors */
-
-typedef real rvec[DIM];
-
-typedef double dvec[DIM];
-
-typedef real matrix[DIM][DIM];
-
-typedef real tensor[DIM][DIM];
-
-typedef int ivec[DIM];
-
-typedef int imatrix[DIM][DIM];
-
-#ifdef __cplusplus
-
-namespace gmx
-{
-
-/*! \brief
- * C++ class for 3D vectors.
- *
- * \tparam ValueType Type
- *
- * This class provides a C++ version of rvec/dvec/ivec that can be put into STL
- * containers etc. It is more or less a drop-in replacement for `rvec` and
- * friends: it can be used in most contexts that accept the equivalent C type.
- * However, there are two cases where explicit conversion is necessary:
- * - An array of these objects needs to be converted with as_vec_array() (or
- * convenience methods like as_rvec_array()).
- * - Passing an RVec as a `const rvec &` parameter to a function needs an
- * explicit call to as_vec(). The implicit conversion should work for this
- * as well, but cppcheck parses the necessary implicit conversion operator
- * incorrectly and MSVC fails to compile code that relies on the implicit
- * conversion, so the explicit method is necessary.
- *
- * For the array conversion to work, the compiler should not add any extra
- * alignment/padding in the layout of this class; that this actually works as
- * intended is tested in the unit tests.
- *
- * \inpublicapi
- */
-template <typename ValueType>
-class BasicVector
-{
- public:
- //! Underlying raw C array type (rvec/dvec/ivec).
- typedef ValueType RawArray[DIM];
-
- //! Constructs default (initialized to 0,0,0) vector.
- BasicVector()
- {
- x_[XX] = 0;
- x_[YY] = 0;
- x_[ZZ] = 0;
- }
- //! Constructs a vector from given values.
- BasicVector(ValueType x, ValueType y, ValueType z)
- {
- x_[XX] = x;
- x_[YY] = y;
- x_[ZZ] = z;
- }
- /*! \brief
- * Constructs a vector from given values.
- *
- * This constructor is not explicit to support implicit conversions
- * that allow, e.g., calling `std::vector<RVec>:``:push_back()` directly
- * with an `rvec` parameter.
- */
- BasicVector(const RawArray x)
- {
- x_[XX] = x[XX];
- x_[YY] = x[YY];
- x_[ZZ] = x[ZZ];
- }
- //! Indexing operator to make the class work as the raw array.
- ValueType &operator[](int i) { return x_[i]; }
- //! Indexing operator to make the class work as the raw array.
- ValueType operator[](int i) const { return x_[i]; }
- // The conversion functions below could more accurately return
- // RawArray &, but this fails with cppcheck and does not solve the
- // issue with MSVC, so as_vec() should be used instead.
- //! Makes BasicVector usable in contexts where a raw C array is expected.
- operator ValueType *() { return x_; }
- //! Makes BasicVector usable in contexts where a raw C array is expected.
- operator const ValueType *() const { return x_; }
- //! Allow inplace addition for BasicVector
- BasicVector<ValueType> &operator+=(const BasicVector<ValueType> &right)
- {
- for (int i = 0; i < DIM; ++i)
- {
- this->x_[i] += right.x_[i];
- }
- return *this;
- }
- //! Allow inplace substraction for BasicVector
- BasicVector<ValueType> &operator-=(const BasicVector<ValueType> &right)
- {
- for (int i = 0; i < DIM; ++i)
- {
- this->x_[i] -= right.x_[i];
- }
- return *this;
- }
- //! Allow vector addition
- BasicVector<ValueType> operator+(const BasicVector<ValueType> &right) const
- {
- BasicVector<ValueType> temp;
- for (int i = 0; i < DIM; ++i)
- {
- temp.x_[i] = this->x_[i] + right.x_[i];
- }
- return temp;
- }
- //! Allow vector substraction
- BasicVector<ValueType> operator-(const BasicVector<ValueType> &right) const
- {
- BasicVector<ValueType> temp;
- for (int i = 0; i < DIM; ++i)
- {
- temp.x_[i] = this->x_[i] - right.x_[i];
- }
- return temp;
- }
- //! Allow vector scalar multiplication (dot product)
- ValueType dot(const BasicVector<ValueType> &right) const
- {
- ValueType temp = 0;
- for (int i = 0; i < DIM; ++i)
- {
- temp += this->x_[i]*right.x_[i];
- }
- return temp;
- }
-
- //! Allow vector vector multiplication (cross product)
- BasicVector<ValueType> cross(const BasicVector<ValueType> &right) const
- {
- BasicVector<ValueType> temp;
- temp.x_[XX] = this->x_[YY]*right.x_[ZZ]-this->x_[ZZ]*right.x_[YY];
- temp.x_[YY] = this->x_[ZZ]*right.x_[XX]-this->x_[XX]*right.x_[ZZ];
- temp.x_[ZZ] = this->x_[XX]*right.x_[YY]-this->x_[YY]*right.x_[XX];
- return temp;
- }
-
- //! Allow vector scaling (vector by scalar multiply)
- BasicVector<ValueType> scale(const ValueType &right) const
- {
- BasicVector<ValueType> temp;
- for (int i = 0; i < DIM; ++i)
- {
- temp.x_[i] = this->x_[i]*right;
- }
- return temp;
- }
-
- //! Return normalized to unit vector
- BasicVector<ValueType> unitv() const
- {
- ValueType length;
- BasicVector<ValueType> temp;
- length = this->length();
- for (int i = 0; i < DIM; ++i)
- {
- temp.x_[i] = this->x_[i]/length;
- }
- return temp;
- }
-
- //! Length^2 of vector
- ValueType sqlength() const
- {
- ValueType temp = 0;
- for (int i = 0; i < DIM; ++i)
- {
- temp += this->x_[i]*this->x_[i];
- }
- return temp;
- }
-
- //! Length of vector
- ValueType length() const
- {
- ValueType temp = 0;
- for (int i = 0; i < DIM; ++i)
- {
- temp += this->x_[i]*this->x_[i];
- }
- return sqrt(temp);
- }
-
- //! cast to RVec
- BasicVector<real> toRVec() const
- {
- BasicVector<real> temp;
- for (int i = 0; i < DIM; ++i)
- {
- temp[i] = static_cast<real>(this->x_[i]);
- }
- return temp;
- }
-
- //! cast to DVec
- BasicVector<double> toDVec() const
- {
- BasicVector<double> temp;
- for (int i = 0; i < DIM; ++i)
- {
- temp[i] = static_cast<double>(this->x_[i]);
- }
- return temp;
- }
-
- //! Converts to a raw C array where implicit conversion does not work.
- RawArray &as_vec() { return x_; }
- //! Converts to a raw C array where implicit conversion does not work.
- const RawArray &as_vec() const { return x_; }
-
- private:
- RawArray x_;
-};
-
-/*! \brief
- * scale for gmx::BasicVector
- */
-template <typename ValueType, typename BasicVector> static inline
-BasicVector scale(BasicVector v, ValueType s)
-{
- BasicVector tmp;
- tmp = v;
- return tmp.scale(s);
-}
-
-/*! \brief
- * unitv for gmx::BasicVector
- */
-template <typename BasicVector> static inline
-BasicVector unitv(BasicVector v)
-{
- return v.unitv();
-}
-
-/*! \brief
- * norm for gmx::BasicVector
- */
-template <typename ValueType> static inline
-ValueType length(BasicVector<ValueType> v)
-{
- return v.length();
-}
-
-/*! \brief
- * norm2 for gmx::BasicVector
- */
-template <typename ValueType> static inline
-ValueType sqlength(BasicVector<ValueType> v)
-{
- return v.sqlength();
-}
-
-/*! \brief
- * cross product for gmx::BasicVector
- */
-template <typename BasicVector> static inline
-BasicVector cross(BasicVector a, BasicVector b)
-{
- return a.cross(b);
-}
-
-/*! \brief
- * dot product for gmx::BasicVector
- */
-template <typename ValueType> static inline
-ValueType dot(BasicVector<ValueType> a, BasicVector<ValueType> b)
-{
- return a.dot(b);
-}
-
-/*! \brief
- * Casts a gmx::BasicVector array into an equivalent raw C array.
- */
-template <typename ValueType> static inline
-typename BasicVector<ValueType>::RawArray *
-as_vec_array(BasicVector<ValueType> *x)
-{
- return reinterpret_cast<typename BasicVector<ValueType>::RawArray *>(x);
-}
-
-/*! \brief
- * Casts a gmx::BasicVector array into an equivalent raw C array.
- */
-template <typename ValueType> static inline
-const typename BasicVector<ValueType>::RawArray *
-as_vec_array(const BasicVector<ValueType> *x)
-{
- return reinterpret_cast<const typename BasicVector<ValueType>::RawArray *>(x);
-}
-
-//! Shorthand for C++ `rvec`-equivalent type.
-typedef BasicVector<real> RVec;
-//! Shorthand for C++ `dvec`-equivalent type.
-typedef BasicVector<double> DVec;
-//! Shorthand for C++ `ivec`-equivalent type.
-typedef BasicVector<int> IVec;
-//! Casts a gmx::RVec array into an `rvec` array.
-static inline rvec *as_rvec_array(RVec *x)
-{
- return as_vec_array(x);
-}
-//! Casts a gmx::DVec array into an `Dvec` array.
-static inline dvec *as_dvec_array(DVec *x)
-{
- return as_vec_array(x);
-}
-//! Casts a gmx::IVec array into an `ivec` array.
-static inline ivec *as_ivec_array(IVec *x)
-{
- return as_vec_array(x);
-}
-
-//! Casts a gmx::RVec array into an `rvec` array.
-static inline const rvec *as_rvec_array(const RVec *x)
-{
- return as_vec_array(x);
-}
-//! Casts a gmx::DVec array into an `dvec` array.
-static inline const dvec *as_dvec_array(const DVec *x)
-{
- return as_vec_array(x);
-}
-//! Casts a gmx::IVec array into an `ivec` array.
-static inline const ivec *as_ivec_array(const IVec *x)
-{
- return as_vec_array(x);
-}
-
-
-} // namespace gmx
-
-#endif
-
-#endif
+++ /dev/null
-/*
- * This file is part of the GROMACS molecular simulation package.
- *
- * Copyright (c) 2012,2013,2014,2015,2016, by the GROMACS development team, led by
- * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
- * and including many others, as listed in the AUTHORS file in the
- * top-level source directory and at http://www.gromacs.org.
- *
- * GROMACS is free software; you can redistribute it and/or
- * modify it under the terms of the GNU Lesser General Public License
- * as published by the Free Software Foundation; either version 2.1
- * of the License, or (at your option) any later version.
- *
- * GROMACS is distributed in the hope that it will be useful,
- * but WITHOUT ANY WARRANTY; without even the implied warranty of
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
- * Lesser General Public License for more details.
- *
- * You should have received a copy of the GNU Lesser General Public
- * License along with GROMACS; if not, see
- * http://www.gnu.org/licenses, or write to the Free Software Foundation,
- * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
- *
- * If you want to redistribute modifications to GROMACS, please
- * consider that scientific software is very special. Version
- * control is crucial - bugs must be traceable. We will be happy to
- * consider code for inclusion in the official distribution, but
- * derived work must not be called official GROMACS. Details are found
- * in the README & COPYING files - if they are missing, get the
- * official version at http://www.gromacs.org.
- *
- * To help us fund GROMACS development, we humbly ask that you cite
- * the research papers on the package. Check out http://www.gromacs.org.
- */
-/*! \libinternal \file
- * \brief
- * Declares OpenMP wrappers to avoid conditional compilation.
- *
- * This module defines wrappers for OpenMP API functions and enables compiling
- * code without conditional compilation even when OpenMP is turned off in the
- * build system.
- * Therefore, OpenMP API functions should always be used through these wrappers
- * and omp.h should never be directly included. Instead, this header should be
- * used whenever OpenMP API functions are needed.
- *
- * \inlibraryapi
- * \ingroup module_utility
- */
-#ifndef GMX_UTILITY_OMP_H
-#define GMX_UTILITY_OMP_H
-
-#include <stdio.h>
-
-#if !GMX_NATIVE_WINDOWS
-/* Ugly hack because the openmp implementation below hacks into the SIMD
- * settings to decide when to use _mm_pause(). This should eventually be
- * changed into proper detection of the intrinsics uses, not SIMD.
- */
-#if GMX_SIMD_X86_SSE2 || GMX_SIMD_X86_SSE4_1 || GMX_SIMD_X86_AVX_128_FMA || \
- GMX_SIMD_X86_AVX_256 || GMX_SIMD_X86_AVX2_256
-# include <xmmintrin.h>
-#endif
-#else
-#include <windows.h>
-#endif
-
-#include "gromacs/utility/basedefinitions.h"
-
-#ifdef __cplusplus
-extern "C"
-{
-#endif
-
-/*! \addtogroup module_utility
- * \{
- */
-
-/*! \brief
- * Returns an integer equal to or greater than the number of threads
- * that would be available if a parallel region without num_threads were
- * defined at that point in the code.
- *
- * Acts as a wrapper for omp_get_max_threads().
- */
-int gmx_omp_get_max_threads(void);
-
-/*! \brief
- * Returns the number of processors available when the function is called.
- *
- * Acts as a wrapper around omp_get_num_procs().
- */
-int gmx_omp_get_num_procs(void);
-
-/*! \brief
- * Returns the thread number of the thread executing within its thread team.
- *
- * Acts as a wrapper for omp_get_thread_num().
- */
-int gmx_omp_get_thread_num(void);
-
-/*! \brief
- * Sets the number of threads in subsequent parallel regions, unless overridden
- * by a num_threads clause.
- *
- * Acts as a wrapper for omp_set_num_threads().
- */
-void gmx_omp_set_num_threads(int num_threads);
-
-/*! \brief
- * Check for externally set thread affinity to avoid conflicts with \Gromacs
- * internal setting.
- *
- * \param[out] message Receives the message to be shown to the user.
- * \returns `true` if we can set thread affinity ourselves.
- *
- * While GNU OpenMP does not set affinity by default, the Intel OpenMP library
- * does. This conflicts with the internal affinity (especially thread-MPI)
- * setting, results in incorrectly locked threads, and causes dreadful performance.
- *
- * The KMP_AFFINITY environment variable is used by Intel, GOMP_CPU_AFFINITY
- * by the GNU compilers (Intel also honors it well). If any of the variables
- * is set, we should honor it and disable the internal pinning.
- * When using Intel OpenMP, we will disable affinity if the user did not set it
- * manually through one of the aforementioned environment variables.
- *
- * Note that the Intel OpenMP affinity disabling will only take effect if this
- * function is called before the OpenMP library gets initialized, which happens
- * when the first call is made into a compilation unit that contains OpenMP
- * pragmas.
- *
- * If this function returns `false`, the caller is responsible to disable the
- * pinning, show the message from \p *message to the user, and free the memory
- * allocated for \p *message.
- * If the return value is `true`, \p *message is NULL.
- */
-gmx_bool gmx_omp_check_thread_affinity(char **message);
-
-/*! \brief
- * Pause for use in a spin-wait loop.
- */
-static gmx_inline void gmx_pause()
-{
-#ifndef _MSC_VER
- /* Ugly hack because the openmp implementation below hacks into the SIMD
- * settings to decide when to use _mm_pause(). This should eventually be
- * changed into proper detection of the intrinsics uses, not SIMD.
- */
-#if (GMX_SIMD_X86_SSE2 || GMX_SIMD_X86_SSE4_1 || GMX_SIMD_X86_AVX_128_FMA || \
- GMX_SIMD_X86_AVX_256 || GMX_SIMD_X86_AVX2_256) && !defined(__MINGW32__)
- /* Replace with tbb::internal::atomic_backoff when/if we use TBB */
- _mm_pause();
-#elif defined __MIC__
- _mm_delay_32(32);
-#else
- /* No wait for unknown architecture */
-#endif
-#else
- YieldProcessor();
-#endif
-}
-
-/*! \} */
-
-#ifdef __cplusplus
-}
-#endif
-
-#endif
+++ /dev/null
-#ifndef NEWFIT_H
-#define NEWFIT_H
-
-#include "gromacs/math/vectypes.h"
-#include "vector"
-
-using gmx::RVec;
-
-namespace NewFit {
-
- void center_frame(std::vector< RVec > &cf_frame, std::vector< int > cf_atoms) {
- if (cf_atoms.empty()) {
- cf_atoms.resize(0);
- for (int i = 0; i < cf_frame.size(); i++) {
- cf_atoms.push_back(i);
- }
- }
- RVec cf_temp;
- cf_temp[0] = 0;
- cf_temp[1] = 0;
- cf_temp[2] = 0;
- for (int i = 0; i < cf_atoms.size(); i++) {
- cf_temp += cf_frame[cf_atoms[i]];
- }
- cf_temp[0] /= cf_atoms.size();
- cf_temp[1] /= cf_atoms.size();
- cf_temp[2] /= cf_atoms.size();
- for (int i = 0; i < cf_atoms.size(); i++) {
- cf_frame[cf_atoms[i]] -= cf_temp;
- }
- }
-
- void implement_fit(std::vector< RVec > &if_frame, RVec if_step, double if_A, double if_B, double if_C) {
- RVec temp;
- for (int i = 0; i < if_frame.size(); i++) {
- if_frame[i] += if_step;
- temp = if_frame[i];
- if_frame[i][0] = cos(if_C) * cos(if_B) * temp[0] + (-sin(if_C)* cos(if_A) + sin(if_B) * cos(if_C) * sin(if_A)) * temp[1] + (sin(if_C) * sin(if_A) + sin(if_B) * cos(if_C) * cos(if_A)) * temp[2];
- if_frame[i][1] = sin(if_C) * cos(if_B) * temp[0] + (cos(if_C) * cos(if_A) + sin(if_C) * sin(if_B) * sin(if_A)) * temp[1] + (-cos(if_C) * sin(if_A) + sin(if_C) * sin(if_B) * cos(if_A)) * temp[2];
- if_frame[i][2] = -sin(if_B) * temp[0] + cos(if_B) * sin(if_A) * temp[1] + cos(if_B) * cos(if_A) * temp[2];
- }
- }
-
- double F(double aix, double aiy, double aiz, double bix, double biy, double biz, double cA, double sA, double cB, double sB, double cC, double sC) {
- return gmx::square(aix + biy*(cA*sC - cC*sA*sB) - biz*(sA*sC + cA*cC*sB) - cB*cC*bix) +
- gmx::square(aiy - biy*(cA*cC + sA*sB*sC) + biz*(cC*sA - cA*sB*sC) - cB*sC*bix) +
- gmx::square(aiz + sB*bix - cA*cB*biz - cB*sA*biy) ;
- }
-
- double Fx(double temp1, double temp2, double temp3, double cB, double sB, double cC, double sC) {
- return 2*sB*temp3 -
- 2*cB*cC*temp1 -
- 2*cB*sC*temp2 ;
- }
-
- double Fy(double temp1, double temp2, double temp3, double cA, double sA, double cB, double sB, double cC, double sC) {
- return 2*(cA*sC - cC*sA*sB)*temp1 -
- 2*(cA*cC + sA*sB*sC)*temp2 -
- 2*cB*sA*temp3 ;
- }
-
- double Fz(double temp1, double temp2, double temp3, double cA, double sA, double cB, double sB, double cC, double sC) {
- return 2*(cC*sA - cA*sB*sC)*temp2 -
- 2*(sA*sC + cA*cC*sB)*temp1 -
- 2*cA*cB*temp3 ;
- }
-
- double FA(double temp1, double temp2, double temp3, double biy_pl_y, double biz_pl_z, double cA, double sA, double cB, double sB, double cC, double sC) {
- return 2*(biy_pl_y*(cC*sA - cA*sB*sC) + biz_pl_z*(cA*cC + sA*sB*sC))*temp2 -
- 2*(cA*cB*biy_pl_y - cB*sA*biz_pl_z)*temp3 -
- 2*(biy_pl_y*(sA*sC + cA*cC*sB) + biz_pl_z*(cA*sC - cC*sA*sB))*temp1 ;
- }
-
- double FB(double temp1, double temp2, double temp3, double bix_pl_x, double biy_pl_y, double biz_pl_z, double cA, double sA, double cB, double sB, double cC, double sC) {
- return 2*(cB*bix_pl_x + sA*sB*biy_pl_y + cA*sB*biz_pl_z)*temp3 -
- 2*(cA*cB*cC*biz_pl_z - cC*sB*bix_pl_x + cB*cC*sA*biy_pl_y)*temp1 -
- 2*(cA*cB*sC*biz_pl_z - sB*sC*bix_pl_x + cB*sA*sC*biy_pl_y)*temp2 ;
- }
-
- double FC(double temp1, double temp2, double bix_pl_x, double biy_pl_y, double biz_pl_z, double cA, double sA, double cB, double sB, double cC, double sC) {
- return 2*(biy_pl_y*(cA*cC + sA*sB*sC) - biz_pl_z*(cC*sA - cA*sB*sC) + cB*sC*bix_pl_x)*temp1 -
- 2*(biz_pl_z*(sA*sC + cA*cC*sB) - biy_pl_y*(cA*sC - cC*sA*sB) + cB*cC*bix_pl_x)*temp2 ;
- }
-
- double dist(std::vector< RVec > alfa, std::vector< RVec > beta, std::vector< std::pair< int, int > > cros) {
- double distance = 0;
- for (int i = 0; i < cros.size(); i++) {
- distance += F(alfa[cros[i].first][0], alfa[cros[i].first][1], alfa[cros[i].first][2], beta[cros[i].second][0], beta[cros[i].second][1], beta[cros[i].second][2], 1, 0, 1, 0, 1, 0);
- }
- return distance;
- }
-
- double dist_old_format(rvec *alfa, rvec *beta, std::vector< std::pair< int, int > > cros) {
- double distance = 0;
- for (int i = 0; i < cros.size(); i++) {
- distance += F(alfa[cros[i].first][0], alfa[cros[i].first][1], alfa[cros[i].first][2], beta[cros[i].second][0], beta[cros[i].second][1], beta[cros[i].second][2], 1, 0, 1, 0, 1, 0);
- }
- return distance;
- }
-
- void new_fit(std::vector< RVec > alfa, std::vector< RVec > &beta, std::vector< std::pair< int, int > > cros, double prec) {
- double bx, by, bz;
- RVec step;
- double x1 = 0, y1 = 0, z1 = 0, A1 = 0, B1 = 0, C1 = 0;
- double x2 = 0, y2 = 0, z2 = 0, A2 = 0, B2 = 0, C2 = 0;
- double cA, sA, cB, sB, cC, sC, temp1, temp2, temp3;
- long double H1 = 0.1, H2 = 0.1, F1 = 0, F2 = 0, H1_temp, H2_temp, H1_const = 2, H2_const = M_PI/90;
- F1 = dist(alfa, beta, cros);
- bool flag = true;
- do {
- F2 = 0;
- for (int i = 0; i < cros.size(); i++) {
- bx = beta[cros[i].second][0] + x1;
- by = beta[cros[i].second][1] + y1;
- bz = beta[cros[i].second][2] + z1;
- cA = cos(A1);
- sA = sin(A1);
- cB = cos(B1);
- sB = sin(B1);
- cC = cos(C1);
- sC = sin(C1);
- temp1 = alfa[cros[i].first][0] + by*(cA*sC - cC*sA*sB) - bz*(sA*sC + cA*cC*sB) - cB*cC*bx;
- temp2 = alfa[cros[i].first][1] - by*(cA*cC + sA*sB*sC) + bz*(cC*sA - cA*sB*sC) - cB*sC*bx;
- temp3 = alfa[cros[i].first][2] + sB*bx - cA*cB*bz - cB*sA*by;
- x2 += Fx(temp1, temp2, temp3, cB, sB, cC, sC);
- y2 += Fy(temp1, temp2, temp3, cA, sA, cB, sB, cC, sC);
- z2 += Fz(temp1, temp2, temp3, cA, sA, cB, sB, cC, sC);
- A2 += FA(temp1, temp2, temp3, by, bz, cA, sA, cB, sB, cC, sC);
- B2 += FB(temp1, temp2, temp3, bx, by, bz, cA, sA, cB, sB, cC, sC);
- C2 += FC(temp1, temp2, bx, by, bz, cA, sA, cB, sB, cC, sC);
- }
- if ((x2 <= prec) && (y2 <= prec) && (z2 <= prec) && (A2 <= prec) && (B2 <= prec) && (C2 <= prec)) {
- flag = false;
- }
- H1_temp = std::max(std::fabs(H1*x2), std::max(std::fabs(H1*y2), std::fabs(H1*z2)));
- if (H1_temp > H1_const) {
- H1 /= 10;
- } else {
- x2 = x1 - H1 * x2;
- y2 = y1 - H1 * y2;
- z2 = z1 - H1 * z2;
- }
- H2_temp = std::max(std::fabs(H2*A2), std::max(std::fabs(H2*B2), std::fabs(H2*C2)));
- if (H2_temp > H2_const) {
- H2 /= 10;
- } else {
- A2 = A1 - H2 * A2;
- B2 = B1 - H2 * B2;
- C2 = C1 - H2 * C2;
- }
- for (int i = 0; i < cros.size(); i++) {
- cA = cos(A2);
- sA = sin(A2);
- cB = cos(B2);
- sB = sin(B2);
- cC = cos(C2);
- sC = sin(C2);
- F2 += F(alfa[cros[i].first][0], alfa[cros[i].first][1], alfa[cros[i].first][2], beta[cros[i].second][0] + x2, beta[cros[i].second][1] + y2, beta[cros[i].second][2] + z2, cA, sA, cB, sB, cC, sC);
- }
- if (F2 <= F1) {
- x1 = x2;
- y1 = y2;
- z1 = z2;
- A1 = A2;
- B1 = B2;
- C1 = C2;
- x2 = 0;
- y2 = 0;
- z2 = 0;
- A2 = 0;
- B2 = 0;
- C2 = 0;
- F1 = F2;
- }
- } while (flag);
- step[0] = x1;
- step[1] = y1;
- step[2] = z1;
- implement_fit(beta, step, A1, B1, C1);
- }
-}
-#endif // NEWFIT_H
#include <iostream>
#include <cstdlib>
-#include <gromacs/utility/gmxomp.h>
-
#include <gromacs/trajectoryanalysis.h>
-#include "gromacs/math/do_fit.h"
-#include "gromacs/utility/smalloc.h"
+#include <gromacs/math/do_fit.h>
+#include <gromacs/utility/smalloc.h>
using namespace gmx;
biod.pnpi.spb.ru / alexxy / gromacs-rcore.git / commitdiff