[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Numerical issues with sin()/cos()
From: |
vlad ionescu |
Subject: |
Numerical issues with sin()/cos() |
Date: |
Mon, 26 Mar 2018 16:01:04 +0000 |
Hello everyone
I am trying to modify fir1.m to create FIRs using the correct sum of
sines or cosines. I ended up making a new .m file, but when I got to
the "%! assert()" lines, I stumbled across some numerical issues. The
.m file was built after a wxMaxima file, which was used to generate
the default results, against which assert() could check, and I noticed
that both sin(), and cos(), seem to not produce (anti-)symmetrical
results. The differences can be more than 1e-15. "format long" is
enabled, but it's the same. Example, lowpass, type I, 33rd order,
wc=0.37, in wxMaxima and Octave (in this order):
0.006248855546642409" "
-0.0151891541716538" "
-0.02000749470712845" "
3.703551396872747*10^-4" "
0.02352639910729612" "
0.01987710435709456" "
-0.0107156702183056" "
-0.03346910762203006" "
-0.01647493965701502" "
0.02756343124616879" "
0.0468059024239567" "
0.006350821104335874" "
-0.06144311077699433" "
-0.07272786646455694" "
0.02972318687964276" "
0.2090466916579446" "
0.3495187814185787" "
0.3495187814185787" "
0.2090466916579446" "
0.02972318687964276" "
-0.07272786646455694" "
-0.06144311077699433" "
0.006350821104335874" "
0.0468059024239567" "
0.02756343124616879" "
-0.01647493965701502" "
-0.03346910762203006" "
-0.0107156702183056" "
0.01987710435709456" "
0.02352639910729612" "
3.703551396872747*10^-4" "
-0.02000749470712845" "
-0.0151891541716538" "
0.006248855546642409
6.24885554664247e-03
-1.51891541716537e-02
-2.00074947071285e-02
3.70355139687233e-04
2.35263991072961e-02
1.98771043570946e-02
-1.07156702183055e-02
-3.34691076220301e-02
-1.64749396570151e-02
2.75634312461687e-02
4.68059024239567e-02
6.35082110433597e-03
-6.14431107769943e-02
-7.27278664645570e-02
2.97231868796425e-02
2.09046691657944e-01
3.49518781418579e-01
3.49518781418577e-01
2.09046691657945e-01
2.97231868796429e-02
-7.27278664645570e-02
-6.14431107769943e-02
6.35082110433577e-03
4.68059024239567e-02
2.75634312461688e-02
-1.64749396570150e-02
-3.34691076220301e-02
-1.07156702183056e-02
1.98771043570945e-02
2.35263991072961e-02
3.70355139687233e-04
-2.00074947071284e-02
-1.51891541716538e-02
6.24885554664247e-03
The simplest generator for this would be (sinc() works the same):
N=33; M=N/2; n=-M:1:M; h=sin(0.37*pi*n)./(pi*n); h(ceil((N+1)/2)=0.37;
Note that the code in Octave directly follows the one in wxMaxima, and
both use the range -M:1:M. There is no cheating in wxMaxima, so
sin(-1)=-sin(1), while in Octave it doesn't seem to be so. What's
more, I am not using big floats in wxMaxima, only default numbers,
float() when necessary.
The assert() that I was aiming for was 1e-15. In the case above, it
fails at h(18):
Location | Observed | Expected | Reason
(18) 0.34952 0.34952 Abs err 1.6653e-15 exceeds tol 1e-15
Is this a known behaviour? Should I just relax the tolerance to 1e-14,
or less? Or am I doing something wrong? I'll attach the .m script as
it is, so far.
Vlad
fir.m
Description: Text Data
- Numerical issues with sin()/cos(),
vlad ionescu <=