Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
4 changes: 3 additions & 1 deletion benchmark/bench_derivative.mojo
Original file line number Diff line number Diff line change
Expand Up @@ -91,7 +91,9 @@ fn poly_func[

fn sin_func[
dtype: DType
](x: Scalar[dtype], args: Optional[List[Scalar[dtype]]]) -> Scalar[dtype] where dtype.is_floating_point():
](x: Scalar[dtype], args: Optional[List[Scalar[dtype]]]) -> Scalar[
dtype
] where dtype.is_floating_point():
return sin(x)


Expand Down
6 changes: 3 additions & 3 deletions pixi.toml
Original file line number Diff line number Diff line change
Expand Up @@ -42,12 +42,12 @@ backend = { name = "pixi-build-mojo", version = "0.*"}
name = "scijo"

[package.host-dependencies]
modular = ">=26.2.0.dev2026022717,<27"
modular = ">=26.2.0.dev2026030605,<27"

[package.build-dependencies]
modular = ">=26.2.0.dev2026022717,<27"
modular = ">=26.2.0.dev2026030605,<27"
numojo = { git = "https://github.com/shivasankarka/NuMojo.git", branch = "pre-0.9"}

[package.run-dependencies]
modular = ">=26.2.0.dev2026022717,<27"
modular = ">=26.2.0.dev2026030605,<27"
numojo = { git = "https://github.com/shivasankarka/NuMojo.git", branch = "pre-0.9"}
5 changes: 0 additions & 5 deletions scijo/fft/__init__.mojo
Original file line number Diff line number Diff line change
Expand Up @@ -11,9 +11,4 @@ The `fft` module provides Fast Fourier Transform operations for complex-valued
arrays. It includes forward and inverse FFT using the Cooley-Tukey algorithm.
"""

from numojo.core.complex import ComplexNDArray, ComplexSIMD
from numojo.core.ndarray import NDArray
from numojo.core.layout import NDArrayShape
from numojo.routines.constants import Constants

from .fastfourier import fft, ifft
67 changes: 47 additions & 20 deletions scijo/fft/fastfourier.mojo
Original file line number Diff line number Diff line change
Expand Up @@ -11,14 +11,18 @@ Forward and inverse Fast Fourier Transform using the Cooley-Tukey radix-2
decimation-in-time algorithm for 1-D complex arrays with power-of-2 lengths.
"""

from math import sin, cos

from numojo.core.complex import ComplexNDArray, ComplexSIMD
from numojo.core.dtype import ComplexDType
from numojo.core.ndarray import NDArray
from numojo.core.layout import NDArrayShape
from numojo.routines.constants import Constants
from numojo.core.indexing import Item

from math import sin, cos
# ===----------------------------------------------------------------------=== #
# FFT
# ===----------------------------------------------------------------------=== #


fn fft[
Expand All @@ -38,25 +42,35 @@ fn fft[
arr: Input complex array to transform. Must be 1-dimensional with length
that is a power of 2.

Raises:
Error: If the input array is not 1-dimensional.
Error: If the array length is not a power of 2.

Returns:
ComplexNDArray containing the FFT of the input array with the same
shape and dtype.

Raises:
Error: If the input array is not 1-dimensional.
Error: If the array length is not a power of 2.
Examples:
```mojo
import numojo as nm
from scijo.fft import fft
from scijo.prelude import *

var arr = nm.linspace[cf32](CScalar[cf32](0, 0), CScalar[cf32](10, 10), num=10)
var fft_arr = fft(arr)
```
"""
if arr.ndim != 1:
raise Error("FFT currently only supports 1D arrays")
raise Error("Scijo [fft]: FFT currently only supports 1D arrays")

var n: Int = arr.shape[0]
if n <= 1:
return ComplexNDArray[dtype](re=arr._re.copy(), im=arr._im.copy())

if (n & (n - 1)) != 0:
raise Error(
"FFT currently only supports arrays with length that is a power"
" of 2"
"Scijo [fft]: FFT currently only supports arrays with length that"
" is a power of 2"
)

var half_size = n // 2
Expand Down Expand Up @@ -109,24 +123,24 @@ fn _ifft_unnormalized[
arr: Input complex array to transform. Must be 1-dimensional with length
that is a power of 2.

Returns:
ComplexNDArray containing the unnormalized inverse FFT of the input array.

Raises:
Error: If the input array is not 1-dimensional.
Error: If the array length is not a power of 2.

Returns:
ComplexNDArray containing the unnormalized inverse FFT of the input array.
"""
if arr.ndim != 1:
raise Error("FFT currently only supports 1D arrays")
raise Error("Scijo [fft]: FFT currently only supports 1D arrays")

var n: Int = arr.shape[0]
if n <= 1:
return ComplexNDArray[dtype](re=arr._re.copy(), im=arr._im.copy())

if (n & (n - 1)) != 0:
raise Error(
"FFT currently only supports arrays with length that is a power"
" of 2"
"Scijo [fft]: FFT currently only supports arrays with length that"
" is a power of 2"
)

var half_size = n // 2
Expand All @@ -144,10 +158,7 @@ fn _ifft_unnormalized[

for k in range(half_size):
var angle = (
2.0
* Constants.pi
* Scalar[dtype.dtype](k)
/ Scalar[dtype.dtype](n)
2.0 * Constants.pi * Scalar[dtype.dtype](k) / Scalar[dtype.dtype](n)
)
var twiddle = ComplexSIMD[dtype](
cos(angle).cast[dtype.dtype](), sin(angle).cast[dtype.dtype]()
Expand All @@ -161,6 +172,11 @@ fn _ifft_unnormalized[
return result^


# ===----------------------------------------------------------------------=== #
# IFFT
# ===----------------------------------------------------------------------=== #


fn ifft[
dtype: ComplexDType = ComplexDType.float64
](arr: ComplexNDArray[dtype]) raises -> ComplexNDArray[
Expand All @@ -178,13 +194,24 @@ fn ifft[
arr: Input complex array to transform. Must be 1-dimensional with length
that is a power of 2.

Raises:
Error: If the input array is not 1-dimensional.
Error: If the array length is not a power of 2.

Returns:
ComplexNDArray containing the IFFT of the input array with the same
shape and dtype.

Raises:
Error: If the input array is not 1-dimensional.
Error: If the array length is not a power of 2.
Examples:
```mojo
import numojo as nm
from scijo.fft import fft, ifft
from scijo.prelude import *

var arr = nm.linspace[cf32](CScalar[cf32](0, 0), CScalar[cf32](10, 10), num=10)
var freq = fft(arr)
var time = ifft(freq)
```
"""
var n: Int = arr.shape[0]
var result = _ifft_unnormalized[dtype](arr)
Expand Down
7 changes: 6 additions & 1 deletion scijo/integrate/__init__.mojo
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,12 @@
The `integrate` module provides tools for numerical integration and quadrature.
It includes adaptive and non-adaptive methods for computing definite integrals,
as well as fixed-sample integration rules for discrete data.

Examples:
```mojo
from scijo.integrate import quad, trapezoid
```
"""

from .quad import quad
from .quadrature import quad
from .fixed_sample import trapezoid, simpson, romb
111 changes: 81 additions & 30 deletions scijo/integrate/fixed_sample.mojo
Original file line number Diff line number Diff line change
Expand Up @@ -12,10 +12,13 @@ Includes the composite trapezoidal rule, Simpson's rule, and Romberg integration
"""

from numojo.core.ndarray import NDArray, NDArrayShape
from numojo.core.error import NumojoError
import numojo as nm


# ===----------------------------------------------------------------------=== #
# Trapezoid
# ===----------------------------------------------------------------------=== #

fn trapezoid[
dtype: DType
](
Expand All @@ -37,37 +40,35 @@ fn trapezoid[
dx: The spacing between sample points. Defaults to 1.0.
axis: The axis along which to integrate. Currently only 1-D is supported.

Raises:
Error: If y is not 1-D.
Error: If y is empty.

Returns:
Definite integral approximated by the trapezoidal rule.
Returns 0.0 for arrays with fewer than 2 elements.

Raises:
Error: If y is not 1-D.
Error: If y is empty.
Examples:
```mojo
import numojo as nm
from scijo.integrate import trapezoid
from scijo.prelude import *

var y = nm.linspace[f64](0.0, 10.0, 100) ** 2 # y = x^2 sampled at 100 points from 0 to 10
var area = trapezoid(y, dx=0.1)
```
"""

if y.ndim != 1:
raise Error(
NumojoError(
category="shape",
message=String(
"Expected y to be 1-D, received ndim={}. Pass a 1-D NDArray"
" for y (e.g. shape (N,))."
).format(y.ndim),
location="trapezoid(y, dx=1.0)",
)
t"Scijo [trapezoid]: Expected y to be 1-D array, received"
t" ndim={{y.ndim}}."
)

if y.size == 0:
raise Error(
NumojoError(
category="value",
message=(
"Cannot integrate over an empty array. Provide a non-empty"
" array for y."
),
location="trapezoid(y, dx=1.0)",
)
t"Scijo [trapezoid]: y.size = 0, Cannot interage over an empty"
t" array."
)

if y.size == 1:
Expand Down Expand Up @@ -103,13 +104,24 @@ fn trapezoid[
x: Array of sample points corresponding to the y values.
axis: The axis along which to integrate. Currently only 1-D is supported.

Raises:
Error: If y or x are not 1-D, or if their sizes differ.
Error: If y is empty.

Returns:
Definite integral approximated by the trapezoidal rule.
Returns 0.0 for arrays with fewer than 2 elements.

Raises:
Error: If y or x are not 1-D, or if their sizes differ.
Error: If y is empty.
Examples:
```mojo
import numojo as nm
from scijo.integrate import trapezoid
from scijo.prelude import *

var x = nm.linspace[f64](0.0, 10.0, 100)
var y = x * x
var area = trapezoid(y, x)
```
"""
if y.ndim != 1:
raise Error(
Expand Down Expand Up @@ -176,6 +188,10 @@ fn trapezoid[
return integral


# ===----------------------------------------------------------------------=== #
# Simpson
# ===----------------------------------------------------------------------=== #

fn simpson[
dtype: DType
](
Expand All @@ -197,11 +213,21 @@ fn simpson[
dx: The spacing between sample points. Defaults to 1.0.
axis: The axis along which to integrate. Currently only 1-D is supported.

Raises:
Error: If y is not 1-D.

Returns:
Definite integral approximated by Simpson's rule.

Raises:
Error: If y is not 1-D.
Examples:
```mojo
import numojo as nm
from scijo.integrate import simpson
from scijo.prelude import *

var y = nm.linspace[f64](0.0, 10.0, 100) ** 2 # y = x^2 sampled at 100 points from 0 to 10
var area = simpson(y, dx=0.1)
```
"""
if y.ndim != 1:
raise Error(
Expand Down Expand Up @@ -251,11 +277,22 @@ fn simpson[
x: Array of sample points corresponding to the y values.
axis: The axis along which to integrate. Currently only 1-D is supported.

Raises:
Error: If y or x are not 1-D, or if their sizes differ.

Returns:
Definite integral approximated by Simpson's rule.

Raises:
Error: If y or x are not 1-D, or if their sizes differ.
Examples:
```mojo
import numojo as nm
from scijo.integrate import simpson
from scijo.prelude import *

var y = nm.linspace[f64](0.0, 10.0, 100) ** 2 # y = x^2 sampled at 100 points from 0 to 10
var x = nm.linspace[f64](0.0, 10.0, 100) # x values corresponding to y
var area = simpson(y, x)
```
"""
if y.ndim != 1:
raise Error(
Expand Down Expand Up @@ -317,6 +354,10 @@ fn simpson[
return integral


# ===----------------------------------------------------------------------=== #
# Romberg
# ===----------------------------------------------------------------------=== #

# TODO: fix the loop implementation.
fn romb[
dtype: DType
Expand All @@ -336,13 +377,23 @@ fn romb[
dx: The spacing between sample points. Defaults to 1.0.
axis: The axis along which to integrate. Currently only 1-D is supported.

Raises:
Error: If y is not 1-D.

Returns:
Definite integral approximated by Romberg integration.

Raises:
Error: If y is not 1-D.
Examples:
```mojo
import numojo as nm
from scijo.integrate import romb
from scijo.prelude import *

var y = nm.linspace[f64](0.0, 10.0, 100) ** 2 # y = x^2 sampled at 100 points from 0 to 10
var area = romb(y, dx=0.1)
```
"""
var maxiter: Int = 10
comptime maxiter: Int = 10
if y.ndim != 1:
raise Error(
NumojoError(
Expand Down
Loading
Loading