diff --git a/CHANGELOG.md b/CHANGELOG.md
index d92591f..61181a1 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -7,7 +7,18 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 ## [Unreleased]
 
+### Changed
+- Errors now surface as normal Python exceptions instead of `PanicException`:
+  - `TypeError` for non-array / non-float32 inputs (with a `.astype(np.float32)` hint).
+  - `ValueError` for wrong input count, rank, or a static dim that doesn't match the model.
+  - `RuntimeError` naming the failing node and op (e.g. `node 'gemm1' (Gemm): ...`) for op
+    and backend failures, including panics raised inside Burn.
+- Rust panic output on stderr is suppressed (the message is in the exception);
+  set `RUST_BACKTRACE=1` to get it back.
+
 ### Fixed
+- Fortran-ordered inputs (e.g. `x.T`) were read in memory order, silently scrambling data.
+  Non-C-contiguous arrays are now copied in logical order.
 - Conv2d / MaxPool2d / AveragePool2d with asymmetric pads no longer average the two
   sides; uneven pads are applied to the input explicitly (exact for `count_include_pad=0`).
 - `auto_pad=SAME_UPPER` / `SAME_LOWER` / `VALID` are now honoured instead of treated as no padding.
diff --git a/src/interpreter/context.rs b/src/interpreter/context.rs
index 4e80143..0f1e2ea 100644
--- a/src/interpreter/context.rs
+++ b/src/interpreter/context.rs
@@ -8,6 +8,9 @@ use onnx_ir::ir::{Argument, ValueSource};
 
 use crate::tensor::{B, FloatPrim, default_device};
 
+/// Ops report failures as a message; the interpreter adds node context.
+pub type OpResult = Result<(), String>;
+
 pub struct ExecutionContext<'w> {
     tensors: HashMap<String, FloatPrim>,
     weights: &'w HashMap<String, FloatPrim>,
@@ -53,4 +56,10 @@ impl<'w> ExecutionContext<'w> {
             ValueSource::Optional => None,
         }
     }
+
+    /// Like `resolve`, but a missing tensor is an error rather than `None`.
+    pub fn require(&self, arg: &Argument) -> Result<FloatPrim, String> {
+        self.resolve(arg)
+            .ok_or_else(|| format!("missing input tensor '{}'", arg.name))
+    }
 }
diff --git a/src/interpreter/mod.rs b/src/interpreter/mod.rs
index 4b31cab..8577f7f 100644
--- a/src/interpreter/mod.rs
+++ b/src/interpreter/mod.rs
@@ -5,18 +5,24 @@ mod context;
 pub mod ops;
 
 use std::collections::HashMap;
+use std::panic::{AssertUnwindSafe, catch_unwind};
 
 use burn_backend::{DType, TensorData, TensorMetadata, backend::ops::FloatTensorOps};
-use numpy::{PyArrayDyn, PyReadonlyArrayDyn};
+use numpy::{
+    PyArrayDyn, PyArrayMethods, PyReadonlyArrayDyn, PyUntypedArray, PyUntypedArrayMethods,
+};
 use onnx_ir::{
     Node, OnnxGraphBuilder,
     batch_norm::{BatchNormConfig, BatchNormalizationNode},
-    ir::{Argument, ValueSource},
+    ir::{ArgType, Argument, ValueSource},
+};
+use pyo3::{
+    exceptions::{PyRuntimeError, PyTypeError, PyValueError},
+    prelude::*,
 };
-use pyo3::prelude::*;
 
 use crate::tensor::{B, FloatPrim, default_device, flex_to_numpy, numpy_to_flex};
-use context::ExecutionContext;
+use context::{ExecutionContext, OpResult};
 
 fn load_weight(arg: &Argument) -> Option<(String, FloatPrim)> {
     if arg.name.is_empty() {
@@ -80,11 +86,95 @@ fn precompute_bn(
     Some((node.name.clone(), scale, offset))
 }
 
+/// A graph input, with whatever shape info the model declares (None = symbolic dim).
+struct InputSpec {
+    name: String,
+    rank: Option<usize>,
+    shape: Option<Vec<Option<usize>>>,
+}
+
+impl InputSpec {
+    fn from_arg(arg: &Argument) -> Self {
+        let (rank, shape) = match &arg.ty {
+            ArgType::Tensor(t) => (Some(t.rank), t.static_shape.clone()),
+            _ => (None, None),
+        };
+        Self {
+            name: arg.name.clone(),
+            rank,
+            shape,
+        }
+    }
+
+    fn check_shape(&self, actual: &[usize]) -> PyResult<()> {
+        let mismatch = match &self.shape {
+            Some(dims) => {
+                dims.len() != actual.len()
+                    || dims
+                        .iter()
+                        .zip(actual)
+                        .any(|(d, a)| d.is_some_and(|d| d != *a))
+            }
+            None => self.rank.is_some_and(|r| r != actual.len()),
+        };
+        if !mismatch {
+            return Ok(());
+        }
+        let expected = match &self.shape {
+            Some(dims) => {
+                let dims: Vec<String> = dims
+                    .iter()
+                    .map(|d| d.map_or("?".to_string(), |d| d.to_string()))
+                    .collect();
+                format!("shape [{}]", dims.join(", "))
+            }
+            None => format!("{} dims", self.rank.unwrap_or_default()),
+        };
+        Err(PyValueError::new_err(format!(
+            "input '{}': expected {expected}, got shape {actual:?}",
+            self.name
+        )))
+    }
+}
+
+/// Accept only float32 ndarrays, with a hint on how to fix anything else.
+fn as_f32_array<'py>(
+    obj: &Bound<'py, PyAny>,
+    name: &str,
+) -> PyResult<PyReadonlyArrayDyn<'py, f32>> {
+    if let Ok(arr) = obj.cast::<PyArrayDyn<f32>>() {
+        return Ok(arr.readonly());
+    }
+    let msg = match obj.cast::<PyUntypedArray>() {
+        Ok(arr) => format!(
+            "input '{name}': expected a float32 array, got {} (use .astype(np.float32))",
+            arr.dtype()
+        ),
+        Err(_) => format!(
+            "input '{name}': expected a numpy array, got {}",
+            obj.get_type().name()?
+        ),
+    };
+    Err(PyTypeError::new_err(msg))
+}
+
+/// Run one node, turning both op errors and backend panics into a message.
+fn run_node(node: &Node, ctx: &mut ExecutionContext) -> OpResult {
+    catch_unwind(AssertUnwindSafe(|| dispatch(node, ctx))).unwrap_or_else(|payload| {
+        let msg = payload
+            .downcast_ref::<String>()
+            .map(String::as_str)
+            .or_else(|| payload.downcast_ref::<&str>().copied())
+            .unwrap_or("unknown panic");
+        Err(msg.to_string())
+    })
+}
+
 #[pyclass]
 pub struct OnnxModel {
     nodes: Vec<Node>,
     weights: HashMap<String, FloatPrim>,
-    input_names: Vec<String>,
+    inputs: Vec<InputSpec>,
     output_names: Vec<String>,
 }
 
@@ -93,40 +183,49 @@ impl OnnxModel {
     fn __call__<'py>(
         &self,
         py: Python<'py>,
-        inputs: Vec<PyReadonlyArrayDyn<'py, f32>>,
+        inputs: Vec<Bound<'py, PyAny>>,
     ) -> PyResult<Vec<Bound<'py, PyArrayDyn<f32>>>> {
-        if inputs.len() != self.input_names.len() {
-            return Err(pyo3::exceptions::PyValueError::new_err(format!(
+        if inputs.len() != self.inputs.len() {
+            return Err(PyValueError::new_err(format!(
                 "expected {} input(s), got {}",
-                self.input_names.len(),
+                self.inputs.len(),
                 inputs.len()
             )));
         }
 
         let mut ctx = ExecutionContext::new(&self.weights);
-        for (name, arr) in self.input_names.iter().zip(inputs) {
-            ctx.insert(name.clone(), numpy_to_flex(&arr));
+        for (spec, obj) in self.inputs.iter().zip(&inputs) {
+            let arr = as_f32_array(obj, &spec.name)?;
+            spec.check_shape(arr.shape())?;
+            ctx.insert(spec.name.clone(), numpy_to_flex(&arr));
         }
 
         for node in &self.nodes {
-            dispatch(node, &mut ctx);
+            run_node(node, &mut ctx).map_err(|msg| {
+                PyRuntimeError::new_err(format!(
+                    "node '{}' ({}): {msg}",
+                    node.name(),
+                    op_type(node)
+                ))
+            })?;
         }
 
         self.output_names
             .iter()
             .map(|name| {
-                let t = ctx
-                    .get(name)
-                    .unwrap_or_else(|| panic!("output '{}' not found", name));
-                Ok(flex_to_numpy(py, t))
+                let t = ctx.get(name).ok_or_else(|| {
+                    PyRuntimeError::new_err(format!("model output '{name}' was never computed"))
+                })?;
+                flex_to_numpy(py, t)
             })
             .collect()
     }
 
     fn __repr__(&self) -> String {
+        let input_names: Vec<&str> = self.inputs.iter().map(|i| i.name.as_str()).collect();
         format!(
             "OnnxModel(inputs={:?}, outputs={:?}, nodes={})",
-            self.input_names,
+            input_names,
             self.output_names,
             self.nodes.len()
         )
@@ -176,7 +275,7 @@ pub fn load_onnx(path: &str) -> Result<OnnxModel, String> {
     }
 
     Ok(OnnxModel {
-        input_names: graph.inputs.iter().map(|a| a.name.clone()).collect(),
+        inputs: graph.inputs.iter().map(InputSpec::from_arg).collect(),
         output_names: graph.outputs.iter().map(|a| a.name.clone()).collect(),
         nodes: graph.nodes,
         weights,
@@ -217,7 +316,7 @@ fn op_type(node: &Node) -> String {
     dbg.split('(').next().unwrap_or(&dbg).to_string()
 }
 
-fn dispatch(node: &Node, ctx: &mut ExecutionContext) {
+fn dispatch(node: &Node, ctx: &mut ExecutionContext) -> OpResult {
     match node {
         Node::Relu(n) => ops::activation::relu(n, ctx),
         Node::Sigmoid(n) => ops::activation::sigmoid(n, ctx),
@@ -239,8 +338,8 @@ fn dispatch(node: &Node, ctx: &mut ExecutionContext) {
         Node::MaxPool2d(n) => ops::pool::max_pool2d(n, ctx),
         Node::AveragePool2d(n) => ops::pool::avg_pool2d(n, ctx),
         Node::GlobalAveragePool(n) => ops::pool::global_avg_pool(n, ctx),
-        Node::Constant(_) => {}
+        Node::Constant(_) => Ok(()),
         // rejected by `is_supported` at load time
-        other => unreachable!("unsupported op '{}'", op_type(other)),
+        other => Err(format!("unsupported op '{}'", op_type(other))),
     }
 }
diff --git a/src/interpreter/ops/activation.rs b/src/interpreter/ops/activation.rs
index b24ba77..7b2175b 100644
--- a/src/interpreter/ops/activation.rs
+++ b/src/interpreter/ops/activation.rs
@@ -1,4 +1,4 @@
-use super::super::context::ExecutionContext;
+use super::super::context::{ExecutionContext, OpResult};
 use crate::tensor::B;
 use burn_backend::backend::ops::{ActivationOps, FloatTensorOps};
 use onnx_ir::{
@@ -6,40 +6,40 @@ use onnx_ir::{
     softmax::SoftmaxNode, tanh::TanhNode,
 };
 
-pub fn relu(node: &ReluNode, ctx: &mut ExecutionContext) {
-    let x = ctx.get(&node.inputs[0].name).expect("relu: missing input");
+pub fn relu(node: &ReluNode, ctx: &mut ExecutionContext) -> OpResult {
+    let x = ctx.require(&node.inputs[0])?;
     ctx.insert(node.outputs[0].name.clone(), B::relu(x));
+    Ok(())
 }
 
-pub fn sigmoid(node: &SigmoidNode, ctx: &mut ExecutionContext) {
-    let x = ctx
-        .get(&node.inputs[0].name)
-        .expect("sigmoid: missing input");
+pub fn sigmoid(node: &SigmoidNode, ctx: &mut ExecutionContext) -> OpResult {
+    let x = ctx.require(&node.inputs[0])?;
     ctx.insert(node.outputs[0].name.clone(), B::sigmoid(x));
+    Ok(())
 }
 
-pub fn tanh(node: &TanhNode, ctx: &mut ExecutionContext) {
-    let x = ctx.get(&node.inputs[0].name).expect("tanh: missing input");
+pub fn tanh(node: &TanhNode, ctx: &mut ExecutionContext) -> OpResult {
+    let x = ctx.require(&node.inputs[0])?;
     ctx.insert(node.outputs[0].name.clone(), B::float_tanh(x));
+    Ok(())
 }
 
-pub fn gelu(node: &GeluNode, ctx: &mut ExecutionContext) {
-    let x = ctx.get(&node.inputs[0].name).expect("gelu: missing input");
+pub fn gelu(node: &GeluNode, ctx: &mut ExecutionContext) -> OpResult {
+    let x = ctx.require(&node.inputs[0])?;
     ctx.insert(node.outputs[0].name.clone(), B::gelu(x));
+    Ok(())
 }
 
-pub fn softmax(node: &SoftmaxNode, ctx: &mut ExecutionContext) {
-    let x = ctx
-        .get(&node.inputs[0].name)
-        .expect("softmax: missing input");
+pub fn softmax(node: &SoftmaxNode, ctx: &mut ExecutionContext) -> OpResult {
+    let x = ctx.require(&node.inputs[0])?;
     let dim = node.config.axis;
     ctx.insert(node.outputs[0].name.clone(), B::softmax(x, dim));
+    Ok(())
 }
 
-pub fn log_softmax(node: &LogSoftmaxNode, ctx: &mut ExecutionContext) {
-    let x = ctx
-        .get(&node.inputs[0].name)
-        .expect("log_softmax: missing input");
+pub fn log_softmax(node: &LogSoftmaxNode, ctx: &mut ExecutionContext) -> OpResult {
+    let x = ctx.require(&node.inputs[0])?;
     let dim = node.config.axis;
     ctx.insert(node.outputs[0].name.clone(), B::log_softmax(x, dim));
+    Ok(())
 }
diff --git a/src/interpreter/ops/conv.rs b/src/interpreter/ops/conv.rs
index 3444a75..e7bb61b 100644
--- a/src/interpreter/ops/conv.rs
+++ b/src/interpreter/ops/conv.rs
@@ -1,4 +1,4 @@
-use super::super::context::ExecutionContext;
+use super::super::context::{ExecutionContext, OpResult};
 use super::padding::{pad2d, resolve_pads, symmetric};
 use crate::tensor::{B, default_device};
 use burn_backend::{
@@ -11,11 +11,9 @@ use onnx_ir::{
     conv2d::Conv2dNode,
 };
 
-pub fn conv2d(node: &Conv2dNode, ctx: &mut ExecutionContext) {
-    let x = ctx.resolve(&node.inputs[0]).expect("conv2d: missing input");
-    let w = ctx
-        .resolve(&node.inputs[1])
-        .expect("conv2d: missing weight");
+pub fn conv2d(node: &Conv2dNode, ctx: &mut ExecutionContext) -> OpResult {
+    let x = ctx.require(&node.inputs[0])?;
+    let w = ctx.require(&node.inputs[1])?;
     let bias = node.inputs.get(2).and_then(|a| ctx.resolve(a));
 
     let cfg = &node.config;
@@ -41,14 +39,13 @@ pub fn conv2d(node: &Conv2dNode, ctx: &mut ExecutionContext) {
     );
 
     ctx.insert(node.outputs[0].name.clone(), B::conv2d(x, w, bias, options));
+    Ok(())
 }
 
 // BatchNorm dispatch — uses pre-fused scale/offset computed at load time.
 // Falls back to full computation if pre-fusion wasn't possible.
-pub fn batch_norm_fused(node: &BatchNormalizationNode, ctx: &mut ExecutionContext) {
-    let x = ctx
-        .resolve(&node.inputs[0])
-        .expect("batch_norm: missing input");
+pub fn batch_norm_fused(node: &BatchNormalizationNode, ctx: &mut ExecutionContext) -> OpResult {
+    let x = ctx.require(&node.inputs[0])?;
 
     let scale_key = format!("{}::scale", node.name);
     let offset_key = format!("{}::offset", node.name);
@@ -60,18 +57,10 @@ pub fn batch_norm_fused(node: &BatchNormalizationNode, ctx: &mut ExecutionContex
         ctx.insert(node.outputs[0].name.clone(), y);
     } else {
         // fallback: full manual BN (slow, but correct)
-        let gamma = ctx
-            .resolve(&node.inputs[1])
-            .expect("batch_norm: missing gamma");
-        let beta = ctx
-            .resolve(&node.inputs[2])
-            .expect("batch_norm: missing beta");
-        let mean = ctx
-            .resolve(&node.inputs[3])
-            .expect("batch_norm: missing mean");
-        let var = ctx
-            .resolve(&node.inputs[4])
-            .expect("batch_norm: missing var");
+        let gamma = ctx.require(&node.inputs[1])?;
+        let beta = ctx.require(&node.inputs[2])?;
+        let mean = ctx.require(&node.inputs[3])?;
+        let var = ctx.require(&node.inputs[4])?;
 
         let eps = match &node.config {
             BatchNormConfig::Static(c) => c.epsilon as f32,
@@ -91,4 +80,5 @@ pub fn batch_norm_fused(node: &BatchNormalizationNode, ctx: &mut ExecutionContex
         let y = B::float_add(y, bcast(beta));
         ctx.insert(node.outputs[0].name.clone(), y);
     }
+    Ok(())
 }
diff --git a/src/interpreter/ops/elementwise.rs b/src/interpreter/ops/elementwise.rs
index 6b536b4..f15cddb 100644
--- a/src/interpreter/ops/elementwise.rs
+++ b/src/interpreter/ops/elementwise.rs
@@ -1,28 +1,32 @@
-use super::super::context::ExecutionContext;
+use super::super::context::{ExecutionContext, OpResult};
 use crate::tensor::B;
 use burn_backend::backend::ops::FloatTensorOps;
 use onnx_ir::arithmetic::{AddNode, DivNode, MulNode, SubNode};
 
-pub fn add(node: &AddNode, ctx: &mut ExecutionContext) {
-    let a = ctx.get(&node.inputs[0].name).expect("add: missing lhs");
-    let b = ctx.get(&node.inputs[1].name).expect("add: missing rhs");
+pub fn add(node: &AddNode, ctx: &mut ExecutionContext) -> OpResult {
+    let a = ctx.require(&node.inputs[0])?;
+    let b = ctx.require(&node.inputs[1])?;
     ctx.insert(node.outputs[0].name.clone(), B::float_add(a, b));
+    Ok(())
 }
 
-pub fn sub(node: &SubNode, ctx: &mut ExecutionContext) {
-    let a = ctx.get(&node.inputs[0].name).expect("sub: missing lhs");
-    let b = ctx.get(&node.inputs[1].name).expect("sub: missing rhs");
+pub fn sub(node: &SubNode, ctx: &mut ExecutionContext) -> OpResult {
+    let a = ctx.require(&node.inputs[0])?;
+    let b = ctx.require(&node.inputs[1])?;
     ctx.insert(node.outputs[0].name.clone(), B::float_sub(a, b));
+    Ok(())
 }
 
-pub fn mul(node: &MulNode, ctx: &mut ExecutionContext) {
-    let a = ctx.get(&node.inputs[0].name).expect("mul: missing lhs");
-    let b = ctx.get(&node.inputs[1].name).expect("mul: missing rhs");
+pub fn mul(node: &MulNode, ctx: &mut ExecutionContext) -> OpResult {
+    let a = ctx.require(&node.inputs[0])?;
+    let b = ctx.require(&node.inputs[1])?;
     ctx.insert(node.outputs[0].name.clone(), B::float_mul(a, b));
+    Ok(())
 }
 
-pub fn div(node: &DivNode, ctx: &mut ExecutionContext) {
-    let a = ctx.get(&node.inputs[0].name).expect("div: missing lhs");
-    let b = ctx.get(&node.inputs[1].name).expect("div: missing rhs");
+pub fn div(node: &DivNode, ctx: &mut ExecutionContext) -> OpResult {
+    let a = ctx.require(&node.inputs[0])?;
+    let b = ctx.require(&node.inputs[1])?;
     ctx.insert(node.outputs[0].name.clone(), B::float_div(a, b));
+    Ok(())
 }
diff --git a/src/interpreter/ops/linear.rs b/src/interpreter/ops/linear.rs
index 58a2a1a..813a642 100644
--- a/src/interpreter/ops/linear.rs
+++ b/src/interpreter/ops/linear.rs
@@ -1,15 +1,11 @@
-use super::super::context::ExecutionContext;
+use super::super::context::{ExecutionContext, OpResult};
 use crate::tensor::{B, default_device};
 use burn_backend::backend::ops::FloatTensorOps;
 use onnx_ir::{gemm::GemmNode, linear::LinearNode};
 
-pub fn linear(node: &LinearNode, ctx: &mut ExecutionContext) {
-    let x = ctx
-        .get(&node.inputs[0].name)
-        .expect("linear: missing input");
-    let w = ctx
-        .get(&node.inputs[1].name)
-        .expect("linear: missing weight");
+pub fn linear(node: &LinearNode, ctx: &mut ExecutionContext) -> OpResult {
+    let x = ctx.require(&node.inputs[0])?;
+    let w = ctx.require(&node.inputs[1])?;
 
     // Gemm layout: w is [out, in], needs transpose before matmul
     // MatMul layout: w is [in, out], use as-is
@@ -43,13 +39,14 @@ pub fn linear(node: &LinearNode, ctx: &mut ExecutionContext) {
     };
 
     ctx.insert(node.outputs[0].name.clone(), y);
+    Ok(())
 }
 
 // General matrix multiply: Y = alpha * A' * B' + beta * C
 // In practice alpha=1, beta=1 for neural net layers.
-pub fn gemm(node: &GemmNode, ctx: &mut ExecutionContext) {
-    let a = ctx.get(&node.inputs[0].name).expect("gemm: missing A");
-    let b = ctx.get(&node.inputs[1].name).expect("gemm: missing B");
+pub fn gemm(node: &GemmNode, ctx: &mut ExecutionContext) -> OpResult {
+    let a = ctx.require(&node.inputs[0])?;
+    let b = ctx.require(&node.inputs[1])?;
 
     let a = if node.config.trans_a != 0 {
         B::float_swap_dims(a, 0, 1)
@@ -95,4 +92,5 @@ pub fn gemm(node: &GemmNode, ctx: &mut ExecutionContext) {
     }
 
     ctx.insert(node.outputs[0].name.clone(), y);
+    Ok(())
 }
diff --git a/src/interpreter/ops/pool.rs b/src/interpreter/ops/pool.rs
index 27f46ad..06a190c 100644
--- a/src/interpreter/ops/pool.rs
+++ b/src/interpreter/ops/pool.rs
@@ -1,4 +1,4 @@
-use super::super::context::ExecutionContext;
+use super::super::context::{ExecutionContext, OpResult};
 use super::padding::{pad2d, resolve_pads, symmetric};
 use crate::tensor::{B, default_device};
 use burn_backend::{
@@ -10,10 +10,8 @@ use onnx_ir::{
     max_pool2d::MaxPool2dNode,
 };
 
-pub fn max_pool2d(node: &MaxPool2dNode, ctx: &mut ExecutionContext) {
-    let x = ctx
-        .resolve(&node.inputs[0])
-        .expect("max_pool2d: missing input");
+pub fn max_pool2d(node: &MaxPool2dNode, ctx: &mut ExecutionContext) -> OpResult {
+    let x = ctx.require(&node.inputs[0])?;
     let cfg = &node.config;
     let pads = resolve_pads(
         &x,
@@ -37,12 +35,11 @@ pub fn max_pool2d(node: &MaxPool2dNode, ctx: &mut ExecutionContext) {
         cfg.ceil_mode,
     );
     ctx.insert(node.outputs[0].name.clone(), y);
+    Ok(())
 }
 
-pub fn avg_pool2d(node: &AveragePool2dNode, ctx: &mut ExecutionContext) {
-    let x = ctx
-        .resolve(&node.inputs[0])
-        .expect("avg_pool2d: missing input");
+pub fn avg_pool2d(node: &AveragePool2dNode, ctx: &mut ExecutionContext) -> OpResult {
+    let x = ctx.require(&node.inputs[0])?;
     let cfg = &node.config;
     let pads = resolve_pads(
         &x,
@@ -80,12 +77,12 @@ pub fn avg_pool2d(node: &AveragePool2dNode, ctx: &mut ExecutionContext) {
         }
     };
     ctx.insert(node.outputs[0].name.clone(), y);
+    Ok(())
 }
 
-pub fn global_avg_pool(node: &GlobalAveragePoolNode, ctx: &mut ExecutionContext) {
-    let x = ctx
-        .resolve(&node.inputs[0])
-        .expect("global_avg_pool: missing input");
+pub fn global_avg_pool(node: &GlobalAveragePoolNode, ctx: &mut ExecutionContext) -> OpResult {
+    let x = ctx.require(&node.inputs[0])?;
     let y = B::adaptive_avg_pool2d(x, [1, 1]);
     ctx.insert(node.outputs[0].name.clone(), y);
+    Ok(())
 }
diff --git a/src/interpreter/ops/reshape.rs b/src/interpreter/ops/reshape.rs
index ae9b0b3..20e2415 100644
--- a/src/interpreter/ops/reshape.rs
+++ b/src/interpreter/ops/reshape.rs
@@ -1,4 +1,4 @@
-use super::super::context::ExecutionContext;
+use super::super::context::{ExecutionContext, OpResult};
 use crate::tensor::B;
 use burn_backend::{TensorMetadata, backend::ops::FloatTensorOps};
 use onnx_ir::{
@@ -7,10 +7,8 @@ use onnx_ir::{
     transpose::TransposeNode,
 };
 
-pub fn reshape(node: &ReshapeNode, ctx: &mut ExecutionContext) {
-    let x = ctx
-        .get(&node.inputs[0].name)
-        .expect("reshape: missing input");
+pub fn reshape(node: &ReshapeNode, ctx: &mut ExecutionContext) -> OpResult {
+    let x = ctx.require(&node.inputs[0])?;
     let orig_shape: Vec<usize> = x.shape().iter().copied().collect();
 
     // shape comes from the attribute (opset < 5) or a constant input (opset 5+)
@@ -18,9 +16,9 @@ pub fn reshape(node: &ReshapeNode, ctx: &mut ExecutionContext) {
         ReshapeInput::Static(shape) => shape.clone(),
         ReshapeInput::Runtime(r) => node.inputs[r.input_index]
             .value()
-            .expect("reshape: shape must be static")
+            .ok_or("shape must be a constant (dynamic shapes are not supported yet)")?
             .to_vec()
-            .expect("reshape: shape must be i64"),
+            .map_err(|e| format!("shape must be int64: {e:?}"))?,
     };
 
     // 0 copies the input dim (allowzero=0); -1 is inferred from what's left
@@ -43,12 +41,11 @@ pub fn reshape(node: &ReshapeNode, ctx: &mut ExecutionContext) {
         node.outputs[0].name.clone(),
         B::float_reshape(x, shape.into()),
     );
+    Ok(())
 }
 
-pub fn flatten(node: &FlattenNode, ctx: &mut ExecutionContext) {
-    let x = ctx
-        .get(&node.inputs[0].name)
-        .expect("flatten: missing input");
+pub fn flatten(node: &FlattenNode, ctx: &mut ExecutionContext) -> OpResult {
+    let x = ctx.require(&node.inputs[0])?;
     let shape: Vec<usize> = x.shape().iter().copied().collect();
     let axis = node.config.axis;
     let outer: usize = shape[..axis].iter().product::<usize>().max(1);
@@ -57,12 +54,11 @@ pub fn flatten(node: &FlattenNode, ctx: &mut ExecutionContext) {
         node.outputs[0].name.clone(),
         B::float_reshape(x, vec![outer, inner].into()),
     );
+    Ok(())
 }
 
-pub fn transpose(node: &TransposeNode, ctx: &mut ExecutionContext) {
-    let x = ctx
-        .get(&node.inputs[0].name)
-        .expect("transpose: missing input");
+pub fn transpose(node: &TransposeNode, ctx: &mut ExecutionContext) -> OpResult {
+    let x = ctx.require(&node.inputs[0])?;
     let rank = x.shape().num_dims();
 
     let perm: Vec<usize> = if node.config.perm.is_empty() {
@@ -72,4 +68,5 @@ pub fn transpose(node: &TransposeNode, ctx: &mut ExecutionContext) {
     };
 
     ctx.insert(node.outputs[0].name.clone(), B::float_permute(x, &perm));
+    Ok(())
 }
diff --git a/src/lib.rs b/src/lib.rs
index 21018c2..a6099ae 100644
--- a/src/lib.rs
+++ b/src/lib.rs
@@ -11,7 +11,7 @@ use interpreter::OnnxModel;
 fn roundtrip<'py>(
     py: Python<'py>,
     arr: PyReadonlyArrayDyn<'py, f32>,
-) -> Bound<'py, PyArrayDyn<f32>> {
+) -> PyResult<Bound<'py, PyArrayDyn<f32>>> {
     let prim = tensor::numpy_to_flex(&arr);
     tensor::flex_to_numpy(py, prim)
 }
@@ -24,6 +24,11 @@ fn load_onnx(path: &str) -> PyResult<OnnxModel> {
 
 #[pymodule]
 fn _burn_python(m: &Bound<'_, PyModule>) -> PyResult<()> {
+    // Panics inside ops are caught and re-raised as Python exceptions that carry the
+    // message, so the default hook's stderr dump is just noise. RUST_BACKTRACE keeps it.
+    if std::env::var_os("RUST_BACKTRACE").is_none() {
+        std::panic::set_hook(Box::new(|_| {}));
+    }
     m.add("__version__", env!("CARGO_PKG_VERSION"))?;
     m.add_function(wrap_pyfunction!(roundtrip, m)?)?;
     m.add_function(wrap_pyfunction!(load_onnx, m)?)?;
diff --git a/src/tensor.rs b/src/tensor.rs
index 71f9836..cdee990 100644
--- a/src/tensor.rs
+++ b/src/tensor.rs
@@ -3,11 +3,11 @@
 // in:  numpy f32 -> TensorData (one copy) -> FlexTensor
 // out: FlexTensor -> TensorData (sync CPU read) -> numpy f32 (one copy)
 
-use burn_backend::{DType, TensorData, backend::ops::FloatTensorOps};
+use burn_backend::{TensorData, backend::ops::FloatTensorOps};
 use burn_flex::{Flex, FlexDevice, FlexTensor};
 use cubecl_common::reader::read_sync;
 use numpy::{PyArray1, PyArrayDyn, PyArrayMethods, PyReadonlyArrayDyn, PyUntypedArrayMethods};
-use pyo3::prelude::*;
+use pyo3::{exceptions::PyRuntimeError, prelude::*};
 
 pub type B = Flex;
 pub type FloatPrim = FlexTensor;
@@ -19,25 +19,28 @@ pub fn default_device() -> FlexDevice {
 // numpy f32 ndarray -> FlexTensor (one copy at the Python boundary)
 pub fn numpy_to_flex(arr: &PyReadonlyArrayDyn<'_, f32>) -> FloatPrim {
     let shape: Vec<usize> = arr.shape().to_vec();
-    let slice = arr
-        .as_slice()
-        .expect("array must be contiguous — call np.ascontiguousarray() first");
-
-    let bytes: Vec<u8> = bytemuck::cast_slice(slice).to_vec();
-    let data = TensorData::from_bytes_vec(bytes, shape, DType::F32);
+    // as_slice() also accepts Fortran order (e.g. x.T), which would be read in memory
+    // order; anything not C-contiguous is gathered in logical order instead
+    let values: Vec<f32> = match arr.as_slice() {
+        Ok(slice) if arr.is_c_contiguous() => slice.to_vec(),
+        _ => arr.as_array().iter().copied().collect(),
+    };
+    let data = TensorData::new(values, shape);
 
     B::float_from_data(data, &default_device())
 }
 
 // FlexTensor -> numpy f32 ndarray (sync read + one copy out)
-pub fn flex_to_numpy<'py>(py: Python<'py>, prim: FloatPrim) -> Bound<'py, PyArrayDyn<f32>> {
+pub fn flex_to_numpy<'py>(
+    py: Python<'py>,
+    prim: FloatPrim,
+) -> PyResult<Bound<'py, PyArrayDyn<f32>>> {
     // for CPU backends the future resolves immediately
-    let data: TensorData = read_sync(B::float_into_data(prim)).expect("float_into_data error");
+    let data: TensorData = read_sync(B::float_into_data(prim))
+        .map_err(|e| PyRuntimeError::new_err(format!("failed to read tensor: {e:?}")))?;
 
     let shape: Vec<usize> = data.shape.iter().copied().collect();
     let floats: &[f32] = bytemuck::cast_slice(data.as_bytes());
 
-    PyArray1::from_slice(py, floats)
-        .reshape(shape)
-        .expect("reshape failed")
+    PyArray1::from_slice(py, floats).reshape(shape)
 }
diff --git a/tests/test_ops.py b/tests/test_ops.py
index 5331674..00d83a0 100644
--- a/tests/test_ops.py
+++ b/tests/test_ops.py
@@ -117,6 +117,43 @@ def reshape(shape):
 except RuntimeError as e:
     check(f"unsupported op raises at load  ({e})", "Exp" in str(e))
 
+def raises(label, fn, exc, *needles):
+    try:
+        fn()
+        check(f"{label}  (no exception)", False)
+    except exc as e:
+        check(f"{label}  ({type(e).__name__}: {e})", all(n in str(e) for n in needles))
+    except BaseException as e:  # e.g. pyo3 PanicException
+        check(f"{label}  (wrong exception {type(e).__name__}: {e})", False)
+
+# ── inference errors ─────────────────────────────────────────────────────────
+
+print("\n=== inference errors ===")
+
+relu = burn.load_onnx(save_model("relu", helper.make_node("Relu", ["x"], ["y"]), ("N", 3)))
+raises("float64 input", lambda: relu([np.zeros((2, 3))]), TypeError, "float32", "float64")
+raises("non-array input", lambda: relu([[1.0, 2.0, 3.0]]), TypeError, "numpy array")
+raises("wrong rank", lambda: relu([np.zeros((2, 3, 1), np.float32)]), ValueError, "[?, 3]")
+raises("wrong static dim", lambda: relu([np.zeros((2, 4), np.float32)]), ValueError, "[?, 3]")
+raises("wrong input count", lambda: relu([]), ValueError, "expected 1 input")
+
+# shape mismatch caught inside the backend: Gemm with an input whose inner dim
+# doesn't match W, declared as fully dynamic so the input check lets it through
+gemm_w = numpy_helper.from_array(np.ones((4, 2), np.float32), "W")
+gemm = burn.load_onnx(save_model(
+    "gemm_dyn", helper.make_node("Gemm", ["x", "W"], ["y"]), ("N", "K"), [gemm_w]))
+raises("backend panic -> RuntimeError with node context",
+       lambda: gemm([np.zeros((2, 3), np.float32)]), RuntimeError, "(Gemm)", "inner dimensions")
+
+# ── non-contiguous inputs ────────────────────────────────────────────────────
+
+print("\n=== non-contiguous inputs ===")
+
+x = np.random.randn(3, 4).astype(np.float32)
+for label, view in [("transposed (F-order)", x.T), ("strided slice", x[:, ::2].T),
+                    ("negative stride", np.ascontiguousarray(x.T)[::-1])]:
+    check(label, np.array_equal(relu([view])[0], np.maximum(view, 0)))
+
 print()
 
 if failures: