EgtNumKernel 1.5a1 : Calcolo degli zeri di polinomi reali e complessi.

This commit is contained in:
Dario Sassi
2014-01-08 14:40:04 +00:00
parent 9dcc85e0db
commit 44b5a1007b
13 changed files with 2334 additions and 0 deletions
+271
View File
@@ -0,0 +1,271 @@
//----------------------------------------------------------------------------
// EgalTech 2013-2013
//----------------------------------------------------------------------------
// File : Complex.cpp Data : 08.01.14 Versione : 1.5a1
// Contenuto : Implementazione classe dei numeri complessi.
//
//
//
// Modifiche : 08.01.14 DS Creazione modulo.
//
//
//----------------------------------------------------------------------------
//--------------------------- Include ----------------------------------------
#include "stdafx.h"
#include "\EgtDev\Include\ENkComplex.h"
//---------------------------- Classe Complex ---------------------------------
Complex
Complex::operator +=( double dVal)
{
this->re += dVal ;
return *this;
}
//----------------------------------------------------------------------------
Complex
Complex::operator +=( Complex& cVal)
{
this->re += cVal.re;
this->im += cVal.im;
return *this;
}
//----------------------------------------------------------------------------
Complex
Complex::operator -=( double dVal)
{
this->re -= dVal ;
return *this ;
}
//----------------------------------------------------------------------------
Complex
Complex::operator -=( Complex& cVal)
{
this->re -= cVal.re ;
this->im -= cVal.im ;
return *this ;
}
//----------------------------------------------------------------------------
Complex
Complex::operator *=( double dVal)
{
this->re *= dVal ;
this->im *= dVal ;
return *this;
}
//----------------------------------------------------------------------------
Complex
Complex::operator /=( double dVal)
{
double dInv ;
dInv = 1.0 / dVal ;
this->re *= dInv ;
this->im *= dInv ;
return *this ;
}
//----------------------------------------------------------------------------
Complex
Complex::operator >>=( int n)
{
this->re = ldexp( this->re, -n) ;
this->im = ldexp( this->im, -n) ;
return *this ;
}
//----------------------------------------------------------------------------
Complex
Complex::operator <<=( int n)
{
this->re = ldexp( this->re, +n) ;
this->im = ldexp( this->im, +n) ;
return *this ;
}
//----------------------------------------------------------------------------
Complex
Complex::operator *=( Complex& cVal)
{
return ( *this = *this * cVal) ;
}
//----------------------------------------------------------------------------
Complex
Complex::operator /=( Complex& cVal)
{
return ( *this *= inv( cVal)) ;
}
//------------------------------ Functions -----------------------------------
// sqrt for Complex
Complex
sqrt( Complex& cVal)
{
Complex z ; // Power 0.5 simple enough
double m ; // to do separate, faster
// than full exp(0.5*log(z))
// just like reals have their
m = mod( cVal) ; // sqrt. Ours gives the one
z.re = sqrt( (m + cVal.re) / 2) ; // with -pi/2 < arg <= +pi/2
z.im = sqrt( (m - cVal.re) / 2) ; // Our log interprets arg
if ( cVal.im < 0.) // as in range -pi to pi,
z.im = - z.im ; // like the atan2 used.
return z ;
}
//----------------------------------------------------------------------------
// log for Complex
Complex
log( Complex& cVal)
{
Complex z ;
z.re = log( m2( cVal)) / 2 ;
z.im = atan2( cVal.im, cVal.re) ;
return z ;
}
//----------------------------------------------------------------------------
Complex
exp( Complex& cVal)
{
Complex ez ;
double m ;
m = exp( cVal.re) ;
ez.re = m * cos( cVal.im) ;
ez.im = m * sin( cVal.im) ;
return ez ;
}
//----------------------------------------------------------------------------
Complex
cosh( Complex& cVal)
{
Complex ez ;
ez = exp( cVal) ;
return ( ( ez + inv(ez)) >> 1) ;
}
//----------------------------------------------------------------------------
Complex
sinh( Complex& cVal)
{
Complex ez ;
ez = exp( cVal) ;
return ( ( ez - inv(ez)) >> 1) ;
}
//----------------------------------------------------------------------------
Complex
tanh( Complex& cVal)
{
Complex e2z ;
e2z = exp( cVal << 1) ;
return ( ( e2z - 1) / ( e2z + 1)) ;
}
//----------------------------------------------------------------------------
Complex
cos( Complex& cVal)
{
return cosh( itimes( cVal)) ;
}
//----------------------------------------------------------------------------
Complex
isin( Complex& cVal)
{
return sinh( itimes( cVal)) ;
}
//----------------------------------------------------------------------------
Complex
sin( Complex& cVal)
{
return -itimes( isin( cVal)) ;
};
//----------------------------------------------------------------------------
Complex
itan( Complex& cVal)
{
return tanh( itimes( cVal)) ;
}
//----------------------------------------------------------------------------
Complex
tan( Complex& cVal)
{
return -itimes( itan( cVal)) ;
}
//----------------------------------------------------------------------------
Complex
acosh( Complex& cVal)
{
return log( cVal + sqrt( cVal * cVal - 1)) ;
}
//----------------------------------------------------------------------------
Complex
asinh( Complex& cVal)
{
return log( cVal + sqrt( cVal * cVal + 1)) ;
}
//----------------------------------------------------------------------------
Complex
atanh( Complex& cVal)
{
return ( log(( 1 + cVal) / ( 1 - cVal)) >> 1) ;
}
//----------------------------------------------------------------------------
Complex
acos( Complex& cVal)
{
return -itimes( acosh( cVal)) ;
}
//----------------------------------------------------------------------------
Complex
asin( Complex& cVal)
{
return -itimes( asinh( itimes( cVal))) ;
}
//----------------------------------------------------------------------------
Complex
atan( Complex& cVal)
{
return -itimes( atanh( itimes( cVal))) ;
}
+20
View File
@@ -0,0 +1,20 @@
//----------------------------------------------------------------------------
// EgalTech 2013-2014
//----------------------------------------------------------------------------
// File : DllMain.h Data : 08.01.14 Versione : 1.5a1
// Contenuto : Prototipi funzioni per uso locale della DLL.
//
//
//
// Modifiche : 08.01.14 DS Creazione modulo.
//
//
//----------------------------------------------------------------------------
#pragma once
#include "/EgtDev/Include/EgtILogger.h"
//-----------------------------------------------------------------------------
ILogger* GetENkLogger( void) ;
+76
View File
@@ -0,0 +1,76 @@
//----------------------------------------------------------------------------
// EgalTech 2013-2013
//----------------------------------------------------------------------------
// File : ENkDllMain.cpp Data : 08.01.14 Versione : 1.5a1
// Contenuto : Inizializzazione della DLL.
//
//
//
// Modifiche : 08.01.14 DS Creazione modulo.
//
//
//----------------------------------------------------------------------------
//--------------------------- Include ----------------------------------------
#include "stdafx.h"
#include "\EgtDev\Include\ENkDllMain.h"
#include "\EgtDev\Include\EgnGetModuleVer.h"
#include "\EgtDev\Include\EgtTrace.h"
//--------------------------- Costanti ----------------------------------------
#if defined( _DEBUG)
const char* ENK_STR = "EgtNumKernelD32.dll ver. " ;
#else
const char* ENK_STR = "EgtNumKernelR32.dll ver. " ;
#endif
const int STR_DIM = 40 ;
//-----------------------------------------------------------------------------
static HINSTANCE s_hModule = NULL ;
static char s_szENkNameVer[STR_DIM] ;
//-----------------------------------------------------------------------------
extern "C" int APIENTRY
DllMain( HMODULE hModule, DWORD dwReason, LPVOID lpReserved)
{
if ( dwReason == DLL_PROCESS_ATTACH) {
s_hModule = hModule ;
EGT_TRACE( "EgtNumKernel.dll Initializing!\n") ;
}
else if ( dwReason == DLL_PROCESS_DETACH) {
s_hModule = NULL ;
EGT_TRACE( "EgtNumKernel.dll Terminating!\n") ;
}
return 1 ;
}
//-----------------------------------------------------------------------------
const char*
GetENkVersion( void)
{
std::string sVer ;
GetModuleVersion( s_hModule, sVer) ;
sprintf_s( s_szENkNameVer, STR_DIM, "%s%s", ENK_STR, sVer.c_str()) ;
return s_szENkNameVer ;
}
//-----------------------------------------------------------------------------
static ILogger* s_pLogger = nullptr ;
//-----------------------------------------------------------------------------
void
SetENkLogger( ILogger* pLogger)
{
s_pLogger = pLogger ;
}
//-----------------------------------------------------------------------------
ILogger*
GetENkLogger( void)
{
return s_pLogger ;
}
BIN
View File
Binary file not shown.
+20
View File
@@ -0,0 +1,20 @@
Microsoft Visual Studio Solution File, Format Version 11.00
# Visual Studio 2010
Project("{8BC9CEB8-8B4A-11D0-8D11-00A0C91BC942}") = "EgtNumKernel", "EgtNumKernel.vcxproj", "{E47BBFD0-36EA-4EC1-9D6D-B05CA9C092C6}"
EndProject
Global
GlobalSection(SolutionConfigurationPlatforms) = preSolution
Debug|Win32 = Debug|Win32
Release|Win32 = Release|Win32
EndGlobalSection
GlobalSection(ProjectConfigurationPlatforms) = postSolution
{E47BBFD0-36EA-4EC1-9D6D-B05CA9C092C6}.Debug|Win32.ActiveCfg = Debug|Win32
{E47BBFD0-36EA-4EC1-9D6D-B05CA9C092C6}.Debug|Win32.Build.0 = Debug|Win32
{E47BBFD0-36EA-4EC1-9D6D-B05CA9C092C6}.Release|Win32.ActiveCfg = Release|Win32
{E47BBFD0-36EA-4EC1-9D6D-B05CA9C092C6}.Release|Win32.Build.0 = Release|Win32
EndGlobalSection
GlobalSection(SolutionProperties) = preSolution
HideSolutionNode = FALSE
EndGlobalSection
EndGlobal
+125
View File
@@ -0,0 +1,125 @@
<?xml version="1.0" encoding="utf-8"?>
<Project DefaultTargets="Build" ToolsVersion="4.0" xmlns="http://schemas.microsoft.com/developer/msbuild/2003">
<ItemGroup Label="ProjectConfigurations">
<ProjectConfiguration Include="Debug|Win32">
<Configuration>Debug</Configuration>
<Platform>Win32</Platform>
</ProjectConfiguration>
<ProjectConfiguration Include="Release|Win32">
<Configuration>Release</Configuration>
<Platform>Win32</Platform>
</ProjectConfiguration>
</ItemGroup>
<PropertyGroup Label="Globals">
<ProjectGuid>{E47BBFD0-36EA-4EC1-9D6D-B05CA9C092C6}</ProjectGuid>
<Keyword>Win32Proj</Keyword>
<RootNamespace>EgtNumKernel</RootNamespace>
</PropertyGroup>
<Import Project="$(VCTargetsPath)\Microsoft.Cpp.Default.props" />
<PropertyGroup Condition="'$(Configuration)|$(Platform)'=='Debug|Win32'" Label="Configuration">
<ConfigurationType>DynamicLibrary</ConfigurationType>
<UseDebugLibraries>true</UseDebugLibraries>
<CharacterSet>Unicode</CharacterSet>
</PropertyGroup>
<PropertyGroup Condition="'$(Configuration)|$(Platform)'=='Release|Win32'" Label="Configuration">
<ConfigurationType>DynamicLibrary</ConfigurationType>
<UseDebugLibraries>false</UseDebugLibraries>
<WholeProgramOptimization>true</WholeProgramOptimization>
<CharacterSet>Unicode</CharacterSet>
</PropertyGroup>
<Import Project="$(VCTargetsPath)\Microsoft.Cpp.props" />
<ImportGroup Label="ExtensionSettings">
</ImportGroup>
<ImportGroup Label="PropertySheets" Condition="'$(Configuration)|$(Platform)'=='Debug|Win32'">
<Import Project="$(UserRootDir)\Microsoft.Cpp.$(Platform).user.props" Condition="exists('$(UserRootDir)\Microsoft.Cpp.$(Platform).user.props')" Label="LocalAppDataPlatform" />
</ImportGroup>
<ImportGroup Label="PropertySheets" Condition="'$(Configuration)|$(Platform)'=='Release|Win32'">
<Import Project="$(UserRootDir)\Microsoft.Cpp.$(Platform).user.props" Condition="exists('$(UserRootDir)\Microsoft.Cpp.$(Platform).user.props')" Label="LocalAppDataPlatform" />
</ImportGroup>
<PropertyGroup Label="UserMacros" />
<PropertyGroup Condition="'$(Configuration)|$(Platform)'=='Debug|Win32'">
<LinkIncremental>true</LinkIncremental>
<OutDir>$(SolutionDir)$(Configuration)$(PlatformArchitecture)\</OutDir>
<IntDir>$(Configuration)$(PlatformArchitecture)\</IntDir>
<TargetName>$(ProjectName)D$(PlatformArchitecture)</TargetName>
</PropertyGroup>
<PropertyGroup Condition="'$(Configuration)|$(Platform)'=='Release|Win32'">
<LinkIncremental>false</LinkIncremental>
<OutDir>$(SolutionDir)$(Configuration)$(PlatformArchitecture)\</OutDir>
<IntDir>$(Configuration)$(PlatformArchitecture)\</IntDir>
<TargetName>$(ProjectName)R$(PlatformArchitecture)</TargetName>
</PropertyGroup>
<ItemDefinitionGroup Condition="'$(Configuration)|$(Platform)'=='Debug|Win32'">
<ClCompile>
<PrecompiledHeader>Use</PrecompiledHeader>
<WarningLevel>Level3</WarningLevel>
<Optimization>Disabled</Optimization>
<PreprocessorDefinitions>WIN32;_DEBUG;_WINDOWS;_USRDLL;I_AM_ENK;%(PreprocessorDefinitions)</PreprocessorDefinitions>
<OpenMPSupport>true</OpenMPSupport>
</ClCompile>
<Link>
<SubSystem>Windows</SubSystem>
<GenerateDebugInformation>true</GenerateDebugInformation>
</Link>
<PostBuildEvent>
<Command>copy $(TargetDir)$(TargetName).pdb \EgtDev\Lib\
copy $(TargetDir)$(TargetName).lib \EgtDev\Lib\
copy $(TargetPath) \EgtProg\Dll</Command>
</PostBuildEvent>
<ResourceCompile>
<PreprocessorDefinitions>_UNICODE;UNICODE;_DEBUG;%(PreprocessorDefinitions)</PreprocessorDefinitions>
</ResourceCompile>
</ItemDefinitionGroup>
<ItemDefinitionGroup Condition="'$(Configuration)|$(Platform)'=='Release|Win32'">
<ClCompile>
<WarningLevel>Level3</WarningLevel>
<PrecompiledHeader>Use</PrecompiledHeader>
<Optimization>MaxSpeed</Optimization>
<FunctionLevelLinking>true</FunctionLevelLinking>
<IntrinsicFunctions>true</IntrinsicFunctions>
<PreprocessorDefinitions>WIN32;NDEBUG;_WINDOWS;_USRDLL;I_AM_ENK;%(PreprocessorDefinitions)</PreprocessorDefinitions>
<OpenMPSupport>true</OpenMPSupport>
</ClCompile>
<Link>
<SubSystem>Windows</SubSystem>
<GenerateDebugInformation>false</GenerateDebugInformation>
<EnableCOMDATFolding>true</EnableCOMDATFolding>
<OptimizeReferences>true</OptimizeReferences>
</Link>
<PostBuildEvent>
<Command>copy $(TargetDir)$(TargetName).pdb \EgtDev\Lib\
copy $(TargetDir)$(TargetName).lib \EgtDev\Lib\
copy $(TargetPath) \EgtProg\Dll</Command>
</PostBuildEvent>
<ResourceCompile>
<PreprocessorDefinitions>_UNICODE;UNICODE;NDEBUG;%(PreprocessorDefinitions)</PreprocessorDefinitions>
</ResourceCompile>
</ItemDefinitionGroup>
<ItemGroup>
<ClCompile Include="Complex.cpp" />
<ClCompile Include="ENkDllMain.cpp" />
<ClCompile Include="JenkinsTraub.cpp" />
<ClCompile Include="PolynomialZeros.cpp" />
<ClCompile Include="stdafx.cpp">
<PrecompiledHeader Condition="'$(Configuration)|$(Platform)'=='Debug|Win32'">Create</PrecompiledHeader>
<PrecompiledHeader Condition="'$(Configuration)|$(Platform)'=='Release|Win32'">Create</PrecompiledHeader>
</ClCompile>
</ItemGroup>
<ItemGroup>
<ClInclude Include="..\Include\EgnGetModuleVer.h" />
<ClInclude Include="..\Include\EgtTargetVer.h" />
<ClInclude Include="..\Include\EgtTrace.h" />
<ClInclude Include="..\Include\ENkDllMain.h" />
<ClInclude Include="..\Include\ENkPolynomialZeros.h" />
<ClInclude Include="DllMain.h" />
<ClInclude Include="JenkinsTraub.h" />
<ClInclude Include="resource.h" />
<ClInclude Include="stdafx.h" />
</ItemGroup>
<ItemGroup>
<ResourceCompile Include="EgtNumKernel.rc" />
</ItemGroup>
<Import Project="$(VCTargetsPath)\Microsoft.Cpp.targets" />
<ImportGroup Label="ExtensionTargets">
</ImportGroup>
</Project>
+68
View File
@@ -0,0 +1,68 @@
<?xml version="1.0" encoding="utf-8"?>
<Project ToolsVersion="4.0" xmlns="http://schemas.microsoft.com/developer/msbuild/2003">
<ItemGroup>
<Filter Include="File di origine">
<UniqueIdentifier>{4FC737F1-C7A5-4376-A066-2A32D752A2FF}</UniqueIdentifier>
<Extensions>cpp;c;cc;cxx;def;odl;idl;hpj;bat;asm;asmx</Extensions>
</Filter>
<Filter Include="File di intestazione">
<UniqueIdentifier>{93995380-89BD-4b04-88EB-625FBE52EBFB}</UniqueIdentifier>
<Extensions>h;hpp;hxx;hm;inl;inc;xsd</Extensions>
</Filter>
<Filter Include="File di risorse">
<UniqueIdentifier>{67DA6AB6-F800-4c08-8B7A-83BB121AAD01}</UniqueIdentifier>
<Extensions>rc;ico;cur;bmp;dlg;rc2;rct;bin;rgs;gif;jpg;jpeg;jpe;resx;tiff;tif;png;wav;mfcribbon-ms</Extensions>
</Filter>
</ItemGroup>
<ItemGroup>
<ClCompile Include="ENkDllMain.cpp">
<Filter>File di origine</Filter>
</ClCompile>
<ClCompile Include="stdafx.cpp">
<Filter>File di origine</Filter>
</ClCompile>
<ClCompile Include="JenkinsTraub.cpp">
<Filter>File di origine</Filter>
</ClCompile>
<ClCompile Include="Complex.cpp">
<Filter>File di origine</Filter>
</ClCompile>
<ClCompile Include="PolynomialZeros.cpp">
<Filter>File di origine</Filter>
</ClCompile>
</ItemGroup>
<ItemGroup>
<ClInclude Include="stdafx.h">
<Filter>File di intestazione</Filter>
</ClInclude>
<ClInclude Include="..\Include\ENkDllMain.h">
<Filter>File di intestazione</Filter>
</ClInclude>
<ClInclude Include="resource.h">
<Filter>File di intestazione</Filter>
</ClInclude>
<ClInclude Include="DllMain.h">
<Filter>File di intestazione</Filter>
</ClInclude>
<ClInclude Include="..\Include\EgnGetModuleVer.h">
<Filter>File di intestazione</Filter>
</ClInclude>
<ClInclude Include="..\Include\EgtTrace.h">
<Filter>File di intestazione</Filter>
</ClInclude>
<ClInclude Include="..\Include\EgtTargetVer.h">
<Filter>File di intestazione</Filter>
</ClInclude>
<ClInclude Include="..\Include\ENkPolynomialZeros.h">
<Filter>File di intestazione</Filter>
</ClInclude>
<ClInclude Include="JenkinsTraub.h">
<Filter>File di intestazione</Filter>
</ClInclude>
</ItemGroup>
<ItemGroup>
<ResourceCompile Include="EgtNumKernel.rc">
<Filter>File di risorse</Filter>
</ResourceCompile>
</ItemGroup>
</Project>
+1390
View File
File diff suppressed because it is too large Load Diff
+90
View File
@@ -0,0 +1,90 @@
//----------------------------------------------------------------------------
// EgalTech 2013-2014
//----------------------------------------------------------------------------
// File : JenkinsTraub.h Data : 08.01.14 Versione : 1.5a1
// Contenuto : Dichiarazione classi per il calcolo degli zeri di polinomi.
//
//
//
// Modifiche : 08.01.14 DS Creazione modulo.
//
//
//----------------------------------------------------------------------------
#pragma once
//----------------------------------------------------------------------------
const int POLY_MAXDEG = 32 ;
//--------------------------- Class Rpoly ------------------------------------
class Rpoly {
public : // methods
int Calculate( const double* op, int degree, double* zeror, double* zeroi) ;
private : // methods
void quad( double a, double b1, double c, double* sr, double* si,
double* lr, double* li) ;
void fxshfr( int l2, int* nz) ;
void quadit( double* uu, double* vv, int* nz) ;
void realit( double sss, int* nz, bool* pIflag) ;
void calcsc( int* type) ;
void nextk( int* type) ;
void newest( int type, double* uu,double* vv) ;
void quadsd( int n, double* u, double* v, double* p, double* q,
double* a, double* b) ;
public : // members
int itercnt ;
private : // members
int n, nn, nmi ;
double sr, si, u, v, a, b, c, d, a1, a2 ;
double a3, a6, a7, e, f, g, h, szr, szi, lzr, lzi ;
double eta, are, mre ;
double p[POLY_MAXDEG+1] ;
double qp[POLY_MAXDEG+1] ;
double k[POLY_MAXDEG+1] ;
double qk[POLY_MAXDEG+1] ;
} ;
//--------------------------- Class Cpoly ------------------------------------
class Cpoly {
public : // methods
int Calculate( const double* opr, const double* opi, int degree, double* zeror, double* zeroi) ;
private : // methods
void noshft( const int l1) ;
void fxshft( const int l2, double* zr, double* zi, bool* pbConv) ;
void vrshft( const int l3, double* zr, double* zi, bool* pbConv) ;
void calct( bool* pbBol) ;
void nexth( bool bBol) ;
void polyev( const int nn, const double sr, const double si, const double pr[], const double pi[],
double qr[], double qi[], double *pvr, double *pvi) ;
double errev( const int nn, const double qr[], const double qi[],
const double ms, const double mp, const double are, const double mre) ;
void cauchy( const int nn, double pt[], double q[], double *fn_val) ;
double scale( const int nn, const double pt[], const double eta, const double infin,
const double smalno, const double base) ;
void cdivid( const double ar, const double ai, const double br, const double bi, double *cr, double *ci) ;
double cmod( const double r, const double i) ;
void mcon( double *eta, double *infiny, double *smalno, double *base) ;
public : // members
int itercnt ;
private : // members
int nn ;
double sr, si, tr, ti, pvr, pvi, are, mre, eta, infin ;
double pr[POLY_MAXDEG+1] ;
double pi[POLY_MAXDEG+1] ;
double hr[POLY_MAXDEG+1] ;
double hi[POLY_MAXDEG+1] ;
double qpr[POLY_MAXDEG+1] ;
double qpi[POLY_MAXDEG+1] ;
double qhr[POLY_MAXDEG+1] ;
double qhi[POLY_MAXDEG+1] ;
double shr[POLY_MAXDEG+1] ;
double shi[POLY_MAXDEG+1] ;
} ;
+236
View File
@@ -0,0 +1,236 @@
//----------------------------------------------------------------------------
// EgalTech 2013-2013
//----------------------------------------------------------------------------
// File : PolynomialZeros.cpp Data : 08.01.14 Versione : 1.5a1
// Contenuto : Funzione per il calcolo degli zeri di polinomi.
//
//
//
// Modifiche : 08.01.14 DS Creazione modulo.
//
//
//----------------------------------------------------------------------------
//--------------------------- Include ----------------------------------------
#include "stdafx.h"
#include "JenkinsTraub.h"
#include "\EgtDev\Include\ENkPolynomialZeros.h"
#include <stdlib.h>
//--------------------------------- Prototipi locali --------------------------------
static void SortRoots( int nNum, double adRoot[]) ;
static void SortRoots( int nNum, Complex acRoot[]) ;
//----------------------------------------------------------------------------
int
PolynomialZeros( int nDegree, double adPoly[], double adRoot[], int* pnIter)
{
int i ;
int j ;
int nZeros ;
double dPreal[POLY_MAXDEG+1] ;
double dZreal[POLY_MAXDEG] ;
double dZcplx[POLY_MAXDEG] ;
Rpoly cRpoly ;
// inizializzo il numero di iterazioni
if ( pnIter != NULL)
*pnIter = 0 ;
// se il coefficiente del grado più alto è zero, diminuisco il grado
while ( nDegree >= 0 && fabs( adPoly[nDegree]) < DBL_EPSILON)
nDegree -- ;
// se il grado è nullo o negativo, errore
if ( nDegree <= 0)
return 0 ;
// verifico di non superare il massimo grado ammesso
if ( nDegree > POLY_MAXDEG)
return 0 ;
// riordino i coefficienti reali
for ( i = 0 ; i <= nDegree ; i++)
dPreal[i] = adPoly[nDegree-i] ;
// calcolo gli zeri
nZeros = cRpoly.Calculate( dPreal, nDegree, dZreal, dZcplx) ;
// assegno gli zeri reali ai parametri di ritorno
for ( i = 0, j = 0 ; i < nZeros ; i++) {
if ( fabs( dZcplx[i]) < 100 * DBL_EPSILON) {
adRoot[j] = dZreal[i] ;
j ++ ;
}
}
nZeros = j ;
// ordino le radici in senso decrescente
SortRoots( nZeros, adRoot) ;
// assegno il numero di iterazioni
if ( pnIter != NULL)
*pnIter = cRpoly.itercnt ;
return nZeros ;
}
//----------------------------------------------------------------------------
int
PolynomialZeros( int nDegree, Complex acPoly[], Complex acRoot[], int* pnIter)
{
bool bCplx ;
int i ;
int nZeros ;
double dPreal[POLY_MAXDEG+1] ;
double dPcplx[POLY_MAXDEG+1] ;
double dZreal[POLY_MAXDEG] ;
double dZcplx[POLY_MAXDEG] ;
Rpoly cRpoly ;
Cpoly cCpoly ;
// inizializzo il numero di iterazioni
if ( pnIter != NULL)
*pnIter = 0 ;
// se il coefficiente del grado più alto è zero, diminuisco il grado
while ( nDegree >= 0 && m2( acPoly[nDegree]) < DBL_EPSILON * DBL_EPSILON)
nDegree -- ;
// se il grado è nullo o negativo, errore
if ( nDegree <= 0)
return 0 ;
// verifico di non superare il massimo grado ammesso
if ( nDegree > POLY_MAXDEG)
return 0 ;
// ricavo i coefficienti reali
for ( i = 0 ; i <= nDegree ; i++)
dPreal[i] = acPoly[nDegree-i].re ;
// ricavo i coefficienti complessi ( e verifico se non nulli)
bCplx = false ;
for ( i = 0 ; i <= nDegree ; i++) {
dPcplx[i] = acPoly[nDegree-i].im ;
if ( fabs( dPcplx[i]) > DBL_EPSILON)
bCplx = true ;
}
// calcolo gli zeri
if ( bCplx)
nZeros = cCpoly.Calculate( dPreal, dPcplx, nDegree, dZreal, dZcplx) ;
else
nZeros = cRpoly.Calculate( dPreal, nDegree, dZreal, dZcplx) ;
// assegno gli zeri ai parametri di ritorno
for ( i = 0 ; i < nZeros ; i++) {
acRoot[i].re = dZreal[i] ;
acRoot[i].im = dZcplx[i] ;
}
// annullo le parti reali e immaginarie molto piccole
for ( i = 0 ; i < nZeros ; i++) {
if ( fabs( acRoot[i].re) < 100 * DBL_EPSILON)
acRoot[i].re = 0 ;
if ( fabs( acRoot[i].im) < 100 * DBL_EPSILON)
acRoot[i].im = 0 ;
}
// ordino le radici in senso decrescente della parte reale
SortRoots( nZeros, acRoot) ;
// assegno il numero di iterazioni
if ( pnIter != NULL)
*pnIter = ( bCplx ? cCpoly.itercnt : cRpoly.itercnt) ;
return nZeros ;
}
//-----------------------------------------------------------------------------
// Confronto tra numeri reali per ordinarli secondo l'ordine crescente
//-----------------------------------------------------------------------------
int
CompareRealRoots( const void* pRoot1, const void* pRoot2)
{
double dRe1 ;
double dRe2 ;
// valori reali
dRe1 = *(double*) pRoot1 ;
dRe2 = *(double*) pRoot2 ;
// se primo maggiore del secondo
if ( dRe1 > dRe2)
return - 1 ;
// se primo minore del secondo
else if ( dRe1 < dRe2)
return + 1 ;
// altrimenti uguali
else
return 0 ;
}
//-----------------------------------------------------------------------------
void
SortRoots( int nNum, double adRoot[])
{
if ( nNum <= 0)
return ;
qsort( adRoot, size_t( nNum), sizeof( double), CompareRealRoots) ;
}
//-----------------------------------------------------------------------------
// Confronto tra numeri complessi per ordinarli secondo l'ordine crescente
// delle parti reali
//-----------------------------------------------------------------------------
int
CompareComplexRoots( const void* pRoot1, const void* pRoot2)
{
double dRe1 ;
double dRe2 ;
double dIm1 ;
double dIm2 ;
// parti reali
dRe1 = Re( *(Complex*) pRoot1) ;
dRe2 = Re( *(Complex*) pRoot2) ;
// se parti reali praticamente uguali
if ( fabs( dRe1 - dRe2) < FLT_MIN) {
// parti immaginarie
dIm1 = Im( *(Complex*) pRoot1) ;
dIm2 = Im( *(Complex*) pRoot2) ;
if ( dIm1 > dIm2)
return - 1 ;
else if ( dIm1 < dIm2)
return + 1 ;
else
return 0 ;
}
// se primo maggiore del secondo
else if ( dRe1 > dRe2)
return - 1 ;
// altrimenti secondo maggiore del primo
else
return + 1 ;
}
//-----------------------------------------------------------------------------
void
SortRoots( int nNum, Complex acRoot[])
{
if ( nNum <= 0)
return ;
qsort( acRoot, size_t( nNum), sizeof( Complex), CompareComplexRoots) ;
}
BIN
View File
Binary file not shown.
+7
View File
@@ -0,0 +1,7 @@
// stdafx.cpp : file di origine che include solo le inclusioni standard
// EgtGeometry.pch sarà l'intestazione precompilata
// stdafx.obj conterrà le informazioni sui tipi precompilati
#include "stdafx.h"
+31
View File
@@ -0,0 +1,31 @@
// stdafx.h : file di inclusione per file di inclusione di sistema standard
// o file di inclusione specifici del progetto utilizzati di frequente, ma
// modificati raramente
//
#pragma once
#include "/EgtDev/Include/EgtTargetVer.h"
#include <stdio.h>
#include <tchar.h>
#include <float.h>
#include <math.h>
// in Debug riconoscimento memory leakage
#if defined( _DEBUG)
#define _CRTDBG_MAP_ALLOC
#include <stdlib.h>
#include <crtdbg.h>
#endif
// in Debug controllo iteratori
#if defined( _DEBUG)
#define _SECURE_SCL 1
#else
#define _SECURE_SCL 0
#endif
#include "/EgtDev/Include/EgtLibVer.h"
#pragma comment(lib, EGTLIBDIR "EgtGeneral" EGTLIBVER ".lib")