//----------------------------------------------------------------------------- // // Author : philippe.billet@noos.fr // // convolution product // // Input: tos-2 x (can be a symbol or expr) // // tos-1 y (can be a symbol or expr) // // // Output: Result on stack // //----------------------------------------------------------------------------- #include "stdafx.h" #include "defs.h" int h1, h2, h3, h4; static void yconvolution(void); void eval_convolution(void) { push(cadr(p1)); eval(); push(caddr(p1)); eval(); convolution(); } void convolution(void) { save(); yconvolution(); restore(); } #define X p1 #define Y p2 static void yconvolution(void) { Y = pop(); X = pop(); push(Y); guess(); p4=pop(); pop(); push(X); guess(); p3=pop(); pop(); if (iszero(X) || iszero(Y)) { push_integer(0); return; } if (car(X) == symbol(DIRAC)) { push(Y); return; } if (car(Y) == symbol(DIRAC)) { push(X); return; } if (car(X) == symbol(DERIVATIVE)) { push_symbol(CONVOLUTION); push(cadr(X)); push(Y); push(caddr(X)); derivative(); list(3); return; } if (car(Y) == symbol(DERIVATIVE)) { push_symbol(CONVOLUTION); push(cadr(Y)); push(X); push(caddr(Y)); derivative(); list(3); return; } if (car(X) == symbol(ADD)) { X = cdr(X); push(car(X)); push(Y); convolution(); X = cdr(X); while (iscons(X)) { push(car(X)); push(Y); convolution(); add(); X = cdr(X); } return; } if (car(Y) == symbol(ADD)) { Y = cdr(Y); push(X); push(car(Y)); convolution(); Y = cdr(Y); while (iscons(Y)) { push(X); push(car(Y)); convolution(); add(); Y = cdr(Y); } return; } if (car(X) == symbol(MULTIPLY) && p3 == p4) { h1 = tos; p5 = cdr(X); while (iscons(p5)) { if (!find(car(p5), p3)) push(car(p5)); p5 = cdr(p5); } h2 = tos; p5 = cdr(X); while (iscons(p5)) { if (find(car(p5), p3)) push(car(p5)); p5 = cdr(p5); } if (tos - h2 == 0) { push(one); push(Y); convolution();} else { multiply_all(tos - h2); X = pop(); push(X); push(Y); convolution(); } multiply_all(tos - h1); return; } if (car(Y) == symbol(MULTIPLY) && p3 == p4) { h3 = tos; p5 = cdr(Y); while (iscons(p5)) { if (!find(car(p5), p3)) push(car(p5)); p5 = cdr(p5); } h4 = tos; p5 = cdr(Y); while (iscons(p5)) { if (find(car(p5), p3)) push(car(p5)); p5 = cdr(p5); } if (tos - h4 == 0) { push(X); push(one); convolution();} else { multiply_all(tos - h4); Y = pop(); push(X); push(Y); convolution(); } multiply_all(tos - h3); return; } // if (exp_flag) { // push(X); // fourier(); // push(Y); // fourier(); // multiply(); // invfourier(); // return; // } push_symbol(CONVOLUTION); push(X); push(Y); list(3); return; } #if SELFTEST static char *s[] = { "convolution(dirac(x),f(x))", "f(x)", "eval(convolution(d(dirac(x),x),f(x)))", "d(f(x),x)", }; void test_convolution(void) { test(__FILE__, s, sizeof s / sizeof (char *)); } #endif