rjhogan / Adept-2

Combined array and automatic differentiation library in C++
http://www.met.reading.ac.uk/clouds/adept/
Apache License 2.0
163 stars 29 forks source link

Wierd behaviour with asin and acos #29

Closed Huzaifg closed 1 year ago

Huzaifg commented 1 year ago

I am seeing weird behavior with asin and acos. Consider the main function in file "summa.cpp"


#include <iostream>
#include <stdint.h>
#include "adept.h"

using adept::adouble;

void simple(std::vector<adouble>& y, std::vector<adouble>& x, adept::Stack& stack)
{

    y[0] = acos(x[2]);
    y[1] = asin(x[2]);
}

int main(int argc, char* argv[]){

    std::vector<double> jac(16);
    adept::Stack stack;
    std::vector<adouble> x(4);
    x[0] = 0;
    x[1] = 0;
    x[2] = 0.5;
    x[3] = 0;

    stack.new_recording();
    std::vector<adouble> y(4);
    simple(y,x, stack);

    stack.independent(&x[0], 4);
    stack.dependent(&y[0], 4);
    stack.jacobian(&jac[0]);

    for(unsigned int i =0; i<16; i++){
        std::cout<<jac[i]<<", ";
    }
    std::cout << std::endl;
    std::cout<<sin(0.5) << std::endl;
    std::cout<<asin(0.5) << std::endl;

}

I would thus expect the following to be printed (computing the jacobian analytically) as Jacobian is stored in column major order

0, 0, 0, 0, 0, 0, 0, 0, -0.479426, 0.877582, 0, 0, 0, 0, 0, 
0.479426
0.479426

However, I get this

0, 0, 0, 0, 0, 0, 0, 0, -1.1547, 1.1547, 0, 0, 0, 0, 0, 0, 
0.479426
0.523599

I use the following MakeFile to compile

CXX = g++
CXXFLAGS = -Wall -g
#Run time seach path for the shared library of adept
LDFLAGS = -Wl,-rpath -Wl,/usr/local/lib
OBJECTS = summa.o
PROGRAM = summa
# Include-file location
INCLUDES = -I/usr/local/include
# Library location and name, plus the math library
LIBS = -L/usr/local/lib -lm -ladept
$(PROGRAM): $(OBJECTS)
    $(CXX) $(CXXFLAGS) $(LDFLAGS) $(OBJECTS) $(LIBS) -o $(PROGRAM)
# Rule to build a normal object file (used to compile all objects in OBJECTS)
%.o: %.cpp
    $(CXX) $(CXXFLAGS) $(LDFLAGS) $(INCLUDES) -c $<

I am running _OS: Arch Linux x8664 Kernel: 6.2.13-arch1-1 g++ (GCC) 13.1.1

I am surely doing something horribly wrong but I am unable to figure out what...

rjhogan commented 1 year ago

This looks correct to me: according to https://en.wikipedia.org/wiki/Inverse_trigonometric_functions, the derivative of asin(x) is 1/sqrt(1-x^2), i.e. 1.1547, which is what Adept gives you. Similarly for acos. [I would recommend not calling your function "sin" if it does something different to the traditional "sin" - this will just confuse users of your code. Adept also supports active and passive arrays, e.g. aVector instead of std::vector, but please read the documentation carefully if you want to do that, particularly when passing aVector objects to and from functions.]

Huzaifg commented 1 year ago

Understood. Thanks alot!