Lezione 11 — Derivazione e integrazione numerica con GNU Octave e GX-TXT

Una curva scientifica può rispondere a domande diverse: quanto vale una grandezza in un determinato istante, quanto rapidamente cambia e quale quantità si accumula nel tempo. Le prime lezioni del corso ci hanno insegnato a costruire, leggere, analizzare e rappresentare Dataset con GNU Octave e GX-TXT. Adesso affrontiamo le ultime due domande attraverso derivazione e integrazione numerica.

Lavoreremo con una curva di concentrazione prodotta da una cinetica del primo ordine: un esperimento sintetico, senza rumore, che permette di confrontare ogni risultato numerico con una soluzione analitica conosciuta. Tutte le elaborazioni scientifiche saranno richieste alle API pubbliche GX-TXT, mentre Octave rimarrà il linguaggio nel quale organizziamo i dati, controlliamo gli errori e scriviamo gli script.

1. Che cosa introduciamo rispetto alle dieci lezioni precedenti

Il corso è progressivo: per comprendere questa lezione non occorre ricominciare dalle variabili o dagli operatori. È però necessario distinguere ciò che già sappiamo dalle nuove operazioni matematiche.

2. L’esperimento: una cinetica di decadimento del primo ordine

Consideriamo una sostanza ipotetica la cui concentrazione diminuisce con velocità proporzionale alla concentrazione stessa:

\[\frac{dC(t)}{dt}=-kC(t)\]

La soluzione analitica dell’equazione è:

\[C(t)=C_0e^{-kt}\]

Scegliamo una concentrazione iniziale 4 mmol/L, una costante 0,055 s⁻¹ e un intervallo di osservazione compreso fra 0 e 60 secondi. Gli esempi rappresentano dati sintetici ideali; non devono essere confusi con misure di laboratorio né con una stima sperimentale di una costante cinetica.

2.1 Perché utilizzare una griglia temporale non uniforme?

Nella Lezione 07 abbiamo generato sinusoidi campionate a intervalli costanti. In questa esercitazione distribuiamo invece i punti in modo non uniforme, concentrandoli all’inizio, dove la concentrazione cambia più rapidamente:

\[t_i=60\left(\frac{i}{100}\right)^{1,4}\mathrm{\,s},\qquad i=0,1,\ldots,100\]

sample_count = 101;
time_s = 60 * (linspace(0, 1, sample_count)' .^ 1.4);
initial_concentration_mmol_L = 4.0;
rate_constant_s_inverse = 0.055;
concentration_mmol_L = initial_concentration_mmol_L * ...
    exp(-rate_constant_s_inverse * time_s);

Le funzioni linspace ed exp, la trasposizione mediante apostrofo e l’operatore elemento per elemento .^ sono già noti. La novità non è la sintassi, ma il significato del vettore time_s: ogni intervallo temporale ha una larghezza potenzialmente differente.

time_steps_s = diff(time_s);
minimum_step_s = min(time_steps_s);
maximum_step_s = max(time_steps_s);
is_strictly_increasing = all(time_steps_s > 0);

Con diff(time_s) ricaviamo gli intervalli fra campioni successivi. Gli estremi del passo temporale sono circa 0,095094 s e 0,838317 s. Si tratta di un’escursione significativa: in un calcolo di derivazione o integrazione sarebbe sbagliato sostituire tutti gli intervalli con un unico passo scelto arbitrariamente.

Apriamo una Session di GX-TXT e l’editor integrato degli script GNU Octave. Gli esempi verranno importati come script di Session con identificativi distinti.

3. La derivata numerica: che cosa significa davvero?

La derivata di una funzione descrive la variazione istantanea rispetto a una coordinata. Se disponiamo di una formula continua, possiamo talvolta derivarla simbolicamente. Se abbiamo soltanto un elenco di campioni, dobbiamo costruire un’approssimazione.

Fra due punti consecutivi, la pendenza della secante è:

\[D_i^{\mathrm{avanti}}=\frac{C_{i+1}-C_i}{t_{i+1}-t_i}\]

Negli intervalli interni, un’altra approssimazione utilizza i due vicini del punto:

\[D_i^{\mathrm{centrata}}\approx\frac{C_{i+1}-C_{i-1}}{t_{i+1}-t_{i-1}}\]

Agli estremi dell’intervallo non è disponibile una coppia di punti su entrambi i lati: si ricorre quindi a secanti unilaterali. È questo il comportamento del kernel numerico utilizzato in questa lezione tramite il servizio GX-TXT.

3.1 Derivata locale contro pendenza di una regressione

La Lezione 10 ha introdotto una retta di regressione, per esempio y = a + b*x. Il coefficiente b è una pendenza stimata sull’insieme delle osservazioni, secondo il modello e il criterio di fitting. Non è, in generale, la derivata locale di una curva non lineare.

La derivata numerica, invece, può restituire una pendenza distinta per ciascun punto. Nel nostro decadimento esponenziale essa è molto negativa all’inizio e si avvicina progressivamente a zero.

3.2 Dimensioni fisiche e segno

Poiché C è espressa in mmol/L e t in secondi, la derivata ha dimensioni mmol/(L·s). Il valore è negativo perché la sostanza diminuisce. Se il fenomeno viene presentato come velocità positiva di scomparsa, occorre definire esplicitamente r(t) = -dC/dt.

Il modello analitico fornisce la verifica indipendente:

\[\frac{dC}{dt}=-kC_0e^{-kt}=-kC(t)\]

dcdt_exact = -rate_constant_s_inverse * concentration_mmol_L;
max_abs_derivative_error = max(abs(dcdt_numeric - dcdt_exact));

4. L’integrale numerico: area e accumulo

L’integrale non descrive una velocità di cambiamento, ma il contributo complessivo di una grandezza lungo un intervallo. Per una concentrazione positiva, l’area sotto la curva è anch’essa positiva.

4.1 Dalla somma alla regola dei trapezi

Con i dati campionati, sostituiamo il tratto di curva fra due punti con un segmento rettilineo. L’area associata è quella di un trapezio:

\[A_i\approx\frac{C_i+C_{i+1}}{2}(t_{i+1}-t_i)\]

Sommando tutti i contributi otteniamo la regola composita dei trapezi:

\[I(T)\approx\sum_{i=0}^{n-2}\frac{C_i+C_{i+1}}{2}(t_{i+1}-t_i)\]

Il fattore (t_(i+1)-t_i) è essenziale. Se gli intervalli sono diversi, usare una sola larghezza per tutti i trapezi introdurrebbe un errore sistematico nella ricostruzione dell’area.

4.2 Integrale definito contro integrale cumulativo

Possiamo voler conoscere soltanto l’integrale sull’intero intervallo, oppure la quantità accumulata fino a ciascun istante:

\[I(T)=\int_{0}^{T}C(t)\,dt\qquad\text{e}\qquad I(t_i)=\int_0^{t_i}C(t)\,dt\]

La media cumulativa della Lezione 09 era una media calcolata su un numero crescente di osservazioni. Un integrale cumulativo è un’altra operazione: somma contributi ponderati dalla larghezza degli intervalli temporali. In generale non coincide né con la somma semplice dei campioni né con una media progressiva.

4.3 Unità di misura dell’integrale

Integrare una concentrazione rispetto al tempo comporta il prodotto delle unità: mmol·s/L. Non è più una concentrazione, né un numero puro. In un contesto applicativo può rappresentare un’esposizione integrata alla concentrazione, ma la sua interpretazione fisica dipende dalla definizione dell’esperimento.

Per la cinetica esponenziale scelta conosciamo anche l’espressione esatta:

\[I(t)=\int_0^t C_0e^{-k\tau}\,d\tau=\frac{C_0}{k}\left(1-e^{-kt}\right)\]

area_exact = initial_concentration_mmol_L / rate_constant_s_inverse * ...
    (1 - exp(-rate_constant_s_inverse * time_s));
exact_total = area_exact(end);

5. L’accuratezza numerica: che cosa stiamo misurando?

Nella Lezione 04 avevamo incontrato l’errore introdotto dall’integrazione temporale di un’equazione differenziale con il metodo di Eulero. Qui incontriamo un errore diverso, pur appartenente alla stessa famiglia di problemi: l’errore di discretizzazione dovuto alla ricostruzione di una derivata o di un’area a partire da punti finiti.

Il confronto con la soluzione analitica permette di misurarlo, senza confonderlo con l’errore di un sensore o con la dispersione delle osservazioni:

\[E_D=\max_i\left|D_{\mathrm{numerico}}(t_i)-D_{\mathrm{esatto}}(t_i)\right|\]

\[E_{I,\mathrm{rel}}=\frac{\left|I_{\mathrm{numerico}}(T)-I_{\mathrm{esatto}}(T)\right|}{\left|I_{\mathrm{esatto}}(T)\right|}\]

La formula dell’errore relativo richiede un riferimento esatto non nullo. Nel nostro esperimento questa condizione è rispettata.

5.1 Perché una griglia più fitta non è sempre la soluzione?

Per dati lisci e privi di rumore, aumentare la densità di campionamento di solito riduce l’errore di approssimazione. Ma i dati sperimentali reali non sono perfetti. Se una misura di concentrazione contiene un errore epsilon, la differenza fra due misure incorpora anche la differenza fra i rispettivi errori:

\[\delta D_i\approx\frac{\varepsilon_{i+1}-\varepsilon_{i-1}}{t_{i+1}-t_{i-1}}\]

Quando l’intervallo al denominatore diminuisce, lo stesso rumore può diventare molto più evidente nella derivata. È un fatto importante in continuità con i segnali rumorosi delle Lezioni 07–09: la derivazione tende a enfatizzare il rumore ad alta frequenza, mentre l’integrazione tende spesso a smorzarne le componenti rapide attraverso l’accumulo.

Filtrare o regolarizzare i dati può essere utile, ma introduce altre scelte e possibili distorsioni. Per non mescolare i fenomeni, questa lezione usa un modello completamente privo di rumore.

6. Le tre nuove operazioni della Numeric API pubblica

NUMERIC_API_v1 non è una nuova libreria privata: è il servizio matematico pubblico GX-TXT già utilizzato nelle lezioni precedenti per operazioni aritmetiche, statistiche e fitting. Nella P008 le operazioni di calcolo differenziale e integrale sono effettivamente richiamabili dall’Octave Script Bridge.

La forma generale della chiamata è:

result = gxtxt_api_call(api_id, operation_name, arguments_cell);

Il primo elemento identifica il servizio, il secondo sceglie l’operazione e il terzo è una cell array contenente gli argomenti posizionali. Le funzioni gxtxt_api_require e gxtxt_api_has_operation permettono di verificare le dipendenze prima del calcolo.

6.1 differentiate: derivare una serie campionata

numeric_api = 'NUMERIC_API_v1';
gxtxt_api_require({numeric_api});
assert(gxtxt_api_has_operation(numeric_api, 'differentiate'));

derivative_result = gxtxt_api_call(numeric_api, 'differentiate', ...
    {time_s, concentration_mmol_L, struct()});
dcdt_numeric = derivative_result.dy_dx(:);

L’API riceve x, y e una struttura di opzioni. I due vettori devono avere la stessa lunghezza, valori numerici finiti e coordinate X strettamente crescenti. La struttura restituita comprende x, dy_dx, n e metadati sul metodo. dy_dx è il vettore delle derivate.

Le funzioni native GNU Octave rimangono utili per studiare la matematica: il kernel effettivo della derivazione usa gradient(y,x). Nello script finale, tuttavia, non eseguiamo un calcolo numerico parallelo e scollegato: chiamiamo il servizio pubblico GX-TXT e analizziamo la sua risposta.

6.2 integrate: integrale definito

integral_result = gxtxt_api_call(numeric_api, 'integrate', ...
    {time_s, concentration_mmol_L, struct()});
area_total = integral_result.value;

L’operazione calcola l’integrale trapezoidale sull’intera griglia e restituisce uno scalare value, insieme a informazioni quali l’intervallo e il metodo. Il kernel Octave corrispondente è trapz(x,y).

6.3 cumulativeIntegrate: un integrale per ogni istante

cumulative_result = gxtxt_api_call(numeric_api, 'cumulativeIntegrate', ...
    {time_s, concentration_mmol_L, struct()});
area_numeric = cumulative_result.integral(:);
assert(abs(area_numeric(end) - area_total) < 1e-8);

Il risultato comprende l’asse x e il vettore integral, allineato alle coordinate originali. Il primo elemento dell’integrale è zero, mentre l’ultimo deve coincidere, entro una tolleranza dichiarata, con l’integrale definito. Il kernel numerico corrispondente in Octave è cumtrapz(x,y).

7. Preparare e pubblicare il Dataset iniziale

Avviamo GX-TXT con una Session nuova. Apriamo Viewers → Octave Scripts, importiamo L11_00_Prepare_Kinetic_Decay.m, assegniamo uno Script ID distinto e salviamo. Dal menu Run → Open Script Runner avviamo il Job.

Lo script genera 101 coppie tempo–concentrazione, le descrive con unità e Variable ID e le pubblica tramite DATASET_PUBLICATION_API_v2/publish_series_batch. La successiva lettura con CANONICAL_DATASET_API_v2/read_range controlla che i valori salvati corrispondano ai valori prodotti.

Il Dataset sorgente ha ID lesson11_decay_observations; la variabile X ha ID time_s, la variabile Y ha ID concentration_mmol_L. Non dobbiamo aprire file di storage né accedere ai registri interni del framework.

8. Un errore scientifico intenzionale: due tempi identici

Il secondo script introduce un difetto nella griglia temporale:

duplicate_time_s = time_s;
duplicate_time_s(32) = duplicate_time_s(31);
time_steps_s = diff(duplicate_time_s);
is_strictly_increasing = all(time_steps_s > 0);
zero_interval_count = sum(time_steps_s == 0);

Il controllo all(time_steps_s > 0) ora restituisce falso, mentre zero_interval_count vale uno. Il problema non è un errore di sintassi; è un dato che viola una condizione matematica richiesta dall’algoritmo.

Apriamo Run → Debug Session script, impostiamo un breakpoint alla riga 35 del file L11_01_Debug_Duplicate_Timestamp.m, dove compare la chiamata all’API, quindi premiamo Start Debug. Nello stato di pausa osserviamo duplicate_time_s, time_steps_s e is_strictly_increasing.

Premendo Continue, la richiesta raggiunge il servizio numerico e produce un errore documentato:

NUMERIC_API_v1/differentiate: math.numeric API v1:
X contains duplicate coordinate values

Questa volta non è uno script dimostrativo a bloccare preventivamente la chiamata: è la Numeric API pubblica che rifiuta la coordinata duplicata.

Non correggiamo automaticamente i dati ordinando il tempo o eliminando campioni. Un riordinamento scorretto potrebbe perdere l’associazione fra coordinate e misure; un’eliminazione potrebbe cancellare informazioni sperimentali valide. La correzione deve nascere dalla verifica della sorgente dei dati.

9. Eseguire la soluzione corretta

Importiamo L11_02_Calculus_Correct.m, salviamolo con un terzo Script ID ed eseguiamolo tramite lo Script Runner. Lo script legge il Dataset originale, esclude valori non finiti e griglie non crescenti, verifica il campionamento non uniforme e invoca le tre operazioni numeriche pubbliche.

Dopo il calcolo, confronta derivate e integrali con le formule analitiche note, applica soglie esplicite di verifica e pubblica quattro nuovi Dataset in un’unica transazione atomica.

9.1 Risultati numerici e interpretazione

La regola dei trapezi restituisce un integrale leggermente maggiore di quello esatto per questa curva esponenziale decrescente, regolare e convessa. Non significa che l’API «introduca una concentrazione in più»: il risultato dipende dall’approssimazione geometrica, dal passo e dalla curvatura della funzione.

10. Visualizzare la derivata e l’integrale in GX-TXT

Il calcolo numerico non termina con una semplice stampa di valori. GX-TXT pubblica quattro Dataset derivati, ciascuno con coordinate, unità, metadati e provenienza. I nomi sono:

Il primo grafico confronta due derivate: entrambe sono negative e si avvicinano a zero con il diminuire della concentrazione. Lo scarto è piccolo rispetto all’ampiezza delle curve, ma resta quantificabile mediante il confronto punto per punto.

Il secondo grafico mostra gli integrali cumulativi: crescono nel tempo, perché la concentrazione rimane positiva. Il valore finale della curva numerica coincide con l’integrale totale restituito da integrate, entro la tolleranza definita.

11. Come lo script collega Octave alle API GX-TXT

Abbiamo incontrato le API di lettura, pubblicazione e grafica nelle lezioni 05–10. Qui le riutilizziamo senza ripeterne tutti i contratti, concentrandoci sulla responsabilità di ciascun servizio.

La sequenza matematica centrale dello script corretto è breve:

derivative_result = gxtxt_api_call(numeric_api, 'differentiate', ...
    {time_s, concentration_mmol_L, struct()});

cumulative_result = gxtxt_api_call(numeric_api, 'cumulativeIntegrate', ...
    {time_s, concentration_mmol_L, struct()});

integral_result = gxtxt_api_call(numeric_api, 'integrate', ...
    {time_s, concentration_mmol_L, struct()});

dcdt_numeric = derivative_result.dy_dx(:);
area_numeric = cumulative_result.integral(:);
area_total = integral_result.value;

I tre risultati hanno forme differenti: il primo fornisce un vettore di derivate, il secondo un vettore di integrali progressivi, il terzo un singolo scalare. Le strutture Octave trasportano anche informazioni sul metodo; non è corretto supporre che qualunque risposta sia direttamente un vettore.

11.1 Perché registriamo un PlotSpec invece di usare plot()?

La funzione grafica nativa di Octave è utile per imparare il linguaggio, ma qui vogliamo una rappresentazione scientifica che resti riconducibile ai risultati pubblicati. Tramite Plot Dataset definiamo esplicitamente assi, unità, Dataset ID e Variable ID. Il grafico rimane quindi una risorsa persistente della Session GX-TXT.

Il metodo già appreso nella Lezione 08 rimane valido: new_spec crea una struttura, validate_spec la verifica, upsert_plot registra la variante, render ne richiede la visualizzazione e artifacts fornisce i riferimenti agli artefatti prodotti.

12. Errori comuni e verifiche essenziali

• La Session non è aperta: lo Script Bridge funziona nel contesto di esecuzione GX-TXT, non come sostituto di un normale terminale Octave.
• Il Dataset iniziale manca: eseguire prima L11_00_Prepare_Kinetic_Decay.m.
• Le coordinate X sono duplicate o disordinate: esaminare diff(time_s); il servizio numerico richiede una griglia valida.
• Il metodo cumulativo non è trovato: attraverso il Bridge usare cumulativeIntegrate con la maiuscola interna.
• Gli ID di uscita sono già presenti: ripetere l’esperimento in una Session nuova; non cancellare direttamente file interni.
• La derivata è molto rumorosa: prima di filtrare, considerare la risoluzione temporale, il rumore di misura e il significato fisico dell’operazione.
• Le curve non compaiono: controllare gli identificativi canonici e validare il PlotSpec prima di renderlo.

13. Esercizi proposti

1. Diminuire il numero dei campioni da 101 a 51. Confrontare l’errore massimo della derivata e quello relativo dell’integrale: quale variazione osservate e perché?
2. Impostare a 1,0 l’esponente usato per costruire il tempo, ottenendo una griglia uniforme. Adeguare il controllo didattico sulla non uniformità e ripetere i calcoli senza cambiare i nomi delle API numeriche.
3. Portare la costante cinetica a 0,10 s⁻¹, aggiornando sia la generazione della curva sia la soluzione analitica di riferimento. Verificare il segno e l’ordine di grandezza della nuova derivata.
4. Ricavare una serie della velocità positiva di scomparsa -dcdt_numeric. Definire un Variable ID e un’unità fisica coerenti prima di pubblicarla.
5. Aggiungere rumore sintetico con una deviazione standard nota. Confrontare la sensibilità della derivata e dell’integrale cumulativo, senza presentare lo scarto numerico come una stima di incertezza sperimentale.
6. Esaminare perché il primo e l’ultimo campione della derivata possono avere un errore diverso da quelli interni. Ripetere il controllo escludendo soltanto dal confronto gli estremi, senza modificare i Dataset originali.

14. Materiale della lezione e riferimenti

I tre script GNU Octave distribuiti con questo articolo sono:

• L11_00_Prepare_Kinetic_Decay.m: genera e pubblica la curva cinetica sintetica.
• L11_01_Debug_Duplicate_Timestamp.m: produce deliberatamente l’errore sulle coordinate per il debugger.
• L11_02_Calculus_Correct.m: esegue il calcolo corretto tramite le API numeriche pubbliche, pubblica i Dataset derivati e crea i grafici persistenti.

Tutti i commenti degli script sono in inglese. Gli esempi sono distribuiti sotto licenza GNU GPL v3 o successiva.

Documentazione pubblica di riferimento:

• Octave Script Bridge API v1
• Numeric API v1
• Canonical Dataset API v2
• Dataset Publication API v2
• Plot Dataset API v4

15. Script GNU Octave completi

Per consultare o copiare gli script senza rendere eccessivamente lunga la pagina, i tre sorgenti sono disponibili nei riquadri espandibili seguenti.

% Lesson 11 - Create a synthetic first-order kinetic decay Dataset.
% Course: GNU Octave with GX-TXT | DEV112 Release 2 + P008
% GitHub release: https://github.com/giuxlab-hub/GX-TXT/releases/tag/v0.9.8.10-dev-112-r2-p008
% P008 source commit (release bundle): f691b0f6e1016cb9af3db30ba3d1c8b5ea91429c
% SPDX-License-Identifier: GPL-3.0-or-later
%
% Run once in a fresh GX-TXT Session with the managed Octave Script Workspace.
% These are synthetic, noise-free observations, not laboratory measurements.

dataset_api = 'CANONICAL_DATASET_API_v2';
publication_api = 'DATASET_PUBLICATION_API_v2';
gxtxt_api_require({dataset_api, publication_api});
assert(gxtxt_api_has_operation(dataset_api, 'list'));
assert(gxtxt_api_has_operation(dataset_api, 'read_range'));
assert(gxtxt_api_has_operation(publication_api, 'publish_series_batch'));

% A nonlinear time transformation creates an ascending, nonuniform grid.
sample_count = 101;
time_s = 60 * (linspace(0, 1, sample_count)' .^ 1.4);
initial_concentration_mmol_L = 4.0;
rate_constant_s_inverse = 0.055;
concentration_mmol_L = initial_concentration_mmol_L * ...
    exp(-rate_constant_s_inverse * time_s);
assert(all(diff(time_s) > 0));
assert(all(isfinite(concentration_mmol_L)));
assert(time_s(1) == 0 && abs(time_s(end)-60) < 1e-12);

source_id = 'lesson11_decay_observations';
present = gxtxt_api_call(dataset_api, 'list', {});
if any(strcmp(present, source_id))
  error('Source Dataset already exists. Create a fresh Session.');
endif

x_variable = struct('variable_id', 'time_s', 'name', 'time_s', ...
  'display_label', 'Elapsed time', 'unit', 's', 'dtype', 'float64', ...
  'value_kind', 'number', 'semantic_role', 'time', ...
  'semantic_source', 'generator');
y_variable = struct('variable_id', 'concentration_mmol_L', ...
  'name', 'concentration_mmol_L', ...
  'display_label', 'Synthetic concentration', 'unit', 'mmol/L', ...
  'dtype', 'float64', 'value_kind', 'number', ...
  'semantic_source', 'generator');

meta = struct();
meta.('simulation.synthetic') = true;
meta.('simulation.kinetic_order') = 1;
meta.('simulation.initial_concentration_mmol_L') = initial_concentration_mmol_L;
meta.('simulation.rate_constant_s_inverse') = rate_constant_s_inverse;
meta.('simulation.nonuniform_sampling') = true;
request = struct('dataset_id', source_id, 'contract', 'table.generic', ...
  'storage', 'auto', 'x_values', time_s, ...
  'y_values', concentration_mmol_L, 'x_variable', x_variable, ...
  'y_variable', y_variable, 'producer_version', '1.0.0', ...
  'source', 'GX-TXT GNU Octave Lesson 11', ...
  'provenance', 'Deterministic synthetic first-order decay', ...
  'notes', 'Nonuniform ascending time grid, noise-free concentrations', ...
  'structured_metadata', meta);
batch = struct('series', {{request}});
receipt = gxtxt_api_call(publication_api, 'publish_series_batch', {batch});
if !receipt.atomic || receipt.output_count != 1
  error('The public publication receipt is not confirmed.');
endif

recovered = gxtxt_api_call(dataset_api, 'read_range', ...
  {source_id, {'time_s', 'concentration_mmol_L'}, 0, sample_count});
assert(max(abs(recovered.time_s.values(:)-time_s)) < 1e-10);
assert(max(abs(recovered.concentration_mmol_L.values(:)-concentration_mmol_L)) < 1e-10);
fprintf('PASS: published %d synthetic kinetic samples as %s.\n', ...
        sample_count, source_id);

% Lesson 11 - Intentional numerical-calculus error for the native debugger.
% Course: GNU Octave with GX-TXT | DEV112 Release 2 + P008
% GitHub release: https://github.com/giuxlab-hub/GX-TXT/releases/tag/v0.9.8.10-dev-112-r2-p008
% P008 source commit (release bundle): f691b0f6e1016cb9af3db30ba3d1c8b5ea91429c
% SPDX-License-Identifier: GPL-3.0-or-later
%
% Run after L11_00. Set a breakpoint immediately before the Numeric API call.
% Inspect duplicate_time_s, time_steps_s, and is_strictly_increasing.
% This script is expected to FAIL inside the public Numeric API, without publishing.

dataset_api = 'CANONICAL_DATASET_API_v2';
numeric_api = 'NUMERIC_API_v1';
gxtxt_api_require({dataset_api, numeric_api});
assert(gxtxt_api_has_operation(dataset_api, 'read_range'));
assert(gxtxt_api_has_operation(numeric_api, 'differentiate'));
observations = gxtxt_api_call(dataset_api, 'read_range', ...
  {'lesson11_decay_observations', {'time_s', 'concentration_mmol_L'}, 0, 101});
time_s = observations.time_s.values(:);
concentration_mmol_L = observations.concentration_mmol_L.values(:);
assert(numel(time_s) == 101 && numel(concentration_mmol_L) == 101);

% INTENTIONAL BUG: one sample repeats the previous timestamp.
% The time order is no longer strictly increasing: public Numeric API must reject it.
duplicate_time_s = time_s;
duplicate_time_s(32) = duplicate_time_s(31);
time_steps_s = diff(duplicate_time_s);
is_strictly_increasing = all(time_steps_s > 0);
zero_interval_count = sum(time_steps_s == 0);
fprintf('Strictly increasing = %d; repeated intervals = %d.\n', ...
        is_strictly_increasing, zero_interval_count);

% BREAKPOINT: inspect duplicate_time_s, time_steps_s and the false flag.
% The Numeric API itself rejects duplicate timestamps; no silent repair.
fprintf('Calling the public GX-TXT derivative with an invalid time grid.\n');
invalid_result = gxtxt_api_call(numeric_api, 'differentiate', ...
    {duplicate_time_s, concentration_mmol_L, struct()});
error('Unexpected success: the Numeric API accepted a duplicate timestamp.');

% Lesson 11 - Numerical differentiation and integration of kinetic data.
% Course: GNU Octave with GX-TXT | DEV112 Release 2 + P008
% GitHub release: https://github.com/giuxlab-hub/GX-TXT/releases/tag/v0.9.8.10-dev-112-r2-p008
% P008 source commit (release bundle): f691b0f6e1016cb9af3db30ba3d1c8b5ea91429c
% SPDX-License-Identifier: GPL-3.0-or-later
%
% Run after L11_00 in the same fresh GX-TXT Session.
% Inputs and outputs are accessed only through public GX-TXT APIs.

dataset_api = 'CANONICAL_DATASET_API_v2';
publication_api = 'DATASET_PUBLICATION_API_v2';
plot_api = 'PLOT_DATASET_API_v4';
numeric_api = 'NUMERIC_API_v1';
gxtxt_api_require({dataset_api, publication_api, plot_api, numeric_api});
required = {
  numeric_api, 'differentiate'; numeric_api, 'integrate';
  numeric_api, 'cumulativeIntegrate';
  dataset_api, 'list'; dataset_api, 'read_range';
  publication_api, 'publish_series_batch';
  plot_api, 'new_spec'; plot_api, 'validate_spec';
  plot_api, 'upsert_plot'; plot_api, 'render'; plot_api, 'artifacts'
};
for k = 1:rows(required)
  if !gxtxt_api_has_operation(required{k,1}, required{k,2})
    error('Missing public operation: %s/%s', required{k,1}, required{k,2});
  endif
endfor

% 1. Read the canonical observations by stable Variable ID.
source_id = 'lesson11_decay_observations';
present = gxtxt_api_call(dataset_api, 'list', {});
if !any(strcmp(present, source_id))
  error('Run L11_00 first in this fresh Session.');
endif
observations = gxtxt_api_call(dataset_api, 'read_range', ...
  {source_id, {'time_s', 'concentration_mmol_L'}, 0, 101});
time_s = observations.time_s.values(:);
concentration_mmol_L = observations.concentration_mmol_L.values(:);
n = numel(time_s);
if n != 101 || numel(concentration_mmol_L) != n || ...
  !all(isfinite(time_s)) || !all(isfinite(concentration_mmol_L)) || ...
  !all(diff(time_s) > 0)
  error('Invalid source data: a finite strictly increasing time grid is required.');
endif

% 2. Inspect the nonuniform grid using standard Octave vector operations.
time_steps_s = diff(time_s);
if max(time_steps_s) - min(time_steps_s) < 1e-8
  error('This experiment requires nonuniform sampling.');
endif

% 3. Use the public GX-TXT numerical CORE, not independent Octave kernels.
% NUMERIC_API_v1/differentiate returns x, dy_dx and method metadata.
derivative_result = gxtxt_api_call(numeric_api, 'differentiate', ...
  {time_s, concentration_mmol_L, struct()});
dcdt_numeric = derivative_result.dy_dx(:);

% The Script Bridge exposes cumulativeIntegrate (camelCase).
% The public Numeric contract describes this as cumulative_integrate.
cumulative_result = gxtxt_api_call(numeric_api, 'cumulativeIntegrate', ...
  {time_s, concentration_mmol_L, struct()});
area_numeric = cumulative_result.integral(:);

% Integrate evaluates the entire range and returns its scalar value.
integral_result = gxtxt_api_call(numeric_api, 'integrate', ...
  {time_s, concentration_mmol_L, struct()});
area_total = integral_result.value;
assert(numel(dcdt_numeric) == n && numel(area_numeric) == n);
assert(all(isfinite(dcdt_numeric)) && all(isfinite(area_numeric)));
assert(isfinite(area_total) && area_total > 0);
assert(max(abs(derivative_result.x(:)-time_s)) < 1e-10);
assert(max(abs(cumulative_result.x(:)-time_s)) < 1e-10);
assert(abs(area_numeric(end) - area_total) < 1e-8);

% 4. Compare with the analytical first-order kinetic reference.
rate_constant_s_inverse = 0.055;
initial_concentration_mmol_L = 4.0;
dcdt_exact = -rate_constant_s_inverse * concentration_mmol_L;
area_exact = initial_concentration_mmol_L / rate_constant_s_inverse * ...
    (1 - exp(-rate_constant_s_inverse * time_s));
exact_total = area_exact(end);
max_abs_derivative_error = max(abs(dcdt_numeric - dcdt_exact));
max_abs_integral_error = max(abs(area_numeric - area_exact));
relative_total_error = abs(area_total - exact_total) / exact_total;
fprintf('Time intervals: minimum %.6f s, maximum %.6f s.\n', ...
        min(diff(time_s)), max(diff(time_s)));
fprintf('Total numerical exposure: %.8f mmol*s/L.\n', area_total);
fprintf('Exact integrated exposure: %.8f mmol*s/L.\n', exact_total);
fprintf('Relative integral error: %.8g.\n', relative_total_error);
fprintf('Max derivative absolute error: %.8g mmol/(L*s).\n', ...
        max_abs_derivative_error);
fprintf('GX-TXT derivative method: %s; integral method: %s.\n', ...
        derivative_result.method, integral_result.method);
if relative_total_error > 0.002 || max_abs_derivative_error > 0.003
  error('Numerical calculus accuracy exceeds lesson tolerance.');
endif

% 5. Publish four scientific series with the true physical dimensions.
derived_ids = {
  'lesson11_derivative_numeric',
  'lesson11_derivative_analytic',
  'lesson11_integral_cumulative',
  'lesson11_integral_analytic'
};
for k = 1:numel(derived_ids)
  if any(strcmp(present, derived_ids{k}))
    error('Output exists: %s. Use a fresh Session.', derived_ids{k});
  endif
endfor
x_variable = struct('variable_id', 'time_s', 'name', 'time_s', ...
  'display_label', 'Elapsed time', 'unit', 's', 'dtype', 'float64', ...
  'value_kind', 'number', 'semantic_role', 'time', ...
  'semantic_source', 'plugin');
rate_variable = struct('variable_id', 'dcdt_mmol_L_s', ...
  'name', 'dcdt_mmol_L_s', 'display_label', 'Concentration derivative', ...
  'unit', 'mmol/(L*s)', 'dtype', 'float64', ...
  'value_kind', 'number', 'semantic_source', 'plugin');
area_variable = struct('variable_id', 'area_mmol_s_L', ...
  'name', 'area_mmol_s_L', 'display_label', 'Cumulative concentration-time integral', ...
  'unit', 'mmol*s/L', 'dtype', 'float64', ...
  'value_kind', 'number', 'semantic_source', 'plugin');
y_values = {dcdt_numeric, dcdt_exact, area_numeric, area_exact};
y_variables = {rate_variable, rate_variable, area_variable, area_variable};
descriptions = {'GX-TXT Numeric API finite-difference derivative', ...
  'Analytical derivative', 'GX-TXT Numeric API cumulative trapezoid', ...
  'Analytical definite integral'};
series = cell(1, numel(derived_ids));
for k = 1:numel(derived_ids)
  metadata = struct();
  metadata.('analysis.method') = descriptions{k};
  metadata.('analysis.source_dataset_id') = source_id;
  metadata.('analysis.nonuniform_grid') = true;
  metadata.('analysis.synthetic_input') = true;
  series{k} = struct('dataset_id', derived_ids{k}, ...
    'contract', 'table.generic', 'storage', 'auto', ...
    'x_values', time_s, 'y_values', y_values{k}, ...
    'x_variable', x_variable, 'y_variable', y_variables{k}, ...
    'producer_version', '1.0.0', ...
    'source', 'GX-TXT GNU Octave Lesson 11', ...
    'source_dataset_ids', {{source_id}}, ...
    'provenance', 'Public GX-TXT NUMERIC_API_v1 calculus on canonical GX-TXT synthetic kinetic data', ...
    'notes', descriptions{k}, 'structured_metadata', metadata);
endfor
batch = struct('series', {series});
receipt = gxtxt_api_call(publication_api, 'publish_series_batch', {batch});
if !receipt.atomic || receipt.output_count != 4
  error('Expected atomic receipt for four derived Datasets.');
endif

% 6. Verify all published values through the canonical public Dataset reader.
for k = 1:numel(derived_ids)
  recovered = gxtxt_api_call(dataset_api, 'read_range', ...
    {derived_ids{k}, {'time_s', y_variables{k}.variable_id}, 0, n});
  saved = recovered.(y_variables{k}.variable_id).values(:);
  assert(numel(saved) == n);
  assert(max(abs(saved - y_values{k})) < 1e-9);
endfor

% 7. Two native, persistent plots preserve Dataset and Variable IDs.
derivative_plot = gxtxt_api_call(plot_api, 'new_spec', ...
  {'plot.dataset.xy', 'Lesson 11: Kinetic derivatives'});
dp = derivative_plot.panels{1};
dp.id = 'derivatives';
dp.title = 'Numerical and analytical concentration derivatives';
dp.axes = struct('x', struct('label', 'Elapsed time [s]', 'scale', 'linear'), ...
  'y_left', struct('label', 'dC/dt [mmol/(L*s)]', 'scale', 'linear'));
dp.series = {
  struct('id', 'numerical', 'label', 'Numerical derivative', ...
    'dataset_id', derived_ids{1}, 'x_variable_id', 'time_s', ...
    'y_variable_id', 'dcdt_mmol_L_s', 'axis', 'left', 'kind', 'line'),
  struct('id', 'analytic', 'label', 'Analytical reference', ...
    'dataset_id', derived_ids{2}, 'x_variable_id', 'time_s', ...
    'y_variable_id', 'dcdt_mmol_L_s', 'axis', 'left', 'kind', 'line')
};
derivative_plot.panels = {dp};
if !gxtxt_api_call(plot_api, 'validate_spec', {derivative_plot})
  error('Numerical derivative PlotSpec invalid.');
endif
derivative_ref = gxtxt_api_call(plot_api, 'upsert_plot', ...
  {'lesson11_kinetic_derivatives', derivative_plot});
gxtxt_api_call(plot_api, 'render', {derivative_ref});
derivative_artifacts = gxtxt_api_call(plot_api, 'artifacts', {derivative_ref});

area_plot = gxtxt_api_call(plot_api, 'new_spec', ...
  {'plot.dataset.xy', 'Lesson 11: Cumulative integral'});
ap = area_plot.panels{1};
ap.id = 'integrals';
ap.title = 'Trapezoidal cumulative area and analytical reference';
ap.axes = struct('x', struct('label', 'Elapsed time [s]', 'scale', 'linear'), ...
  'y_left', struct('label', 'Integral C dt [mmol*s/L]', 'scale', 'linear'));
ap.series = {
  struct('id', 'numerical', 'label', 'Numerical trapezoid', ...
    'dataset_id', derived_ids{3}, 'x_variable_id', 'time_s', ...
    'y_variable_id', 'area_mmol_s_L', 'axis', 'left', 'kind', 'line'),
  struct('id', 'analytic', 'label', 'Analytical reference', ...
    'dataset_id', derived_ids{4}, 'x_variable_id', 'time_s', ...
    'y_variable_id', 'area_mmol_s_L', 'axis', 'left', 'kind', 'line')
};
area_plot.panels = {ap};
if !gxtxt_api_call(plot_api, 'validate_spec', {area_plot})
  error('Numerical integration PlotSpec invalid.');
endif
area_ref = gxtxt_api_call(plot_api, 'upsert_plot', ...
  {'lesson11_kinetic_integral', area_plot});
gxtxt_api_call(plot_api, 'render', {area_ref});
area_artifacts = gxtxt_api_call(plot_api, 'artifacts', {area_ref});
fprintf('PASS: public NUMERIC_API_v1 calculus, 4 derived Datasets, 2 native plots, read-back verified.\n');
disp(derivative_artifacts);
disp(area_artifacts);

Questa esercitazione conclude il lavoro su derivate e integrali numerici e prepara all’impiego dei metodi di calcolo in applicazioni più complesse, nelle quali la scelta della griglia, delle unità e del modello numerico è parte integrante dell’interpretazione scientifica.