forked from animatedread/Warrior_EA
ea9d86b3ee
CaclHiddenGradient computed this layer's gradient as matrix_w[(outputs+1)*i + k] against a buffer whose actual layout (one row per NEXT-layer neuron, stride inputs+1) makes the correct read matrix_w[k*(inputs+1) + i]: the transpose for square layers, and for the non-square boundaries this EA actually builds (tapered stacks, the 3-neuron head) a mis-strided walk that ran past the buffer end - garbage on OpenCL, zeroed reads on the CPU DLL, so the tiers did not even agree with each other. Every gradient crossing a dense boundary on its way down - the entire learning signal reaching the BN/conv/LSTM front ends - passed through a fixed wrong matrix: feedback-alignment dynamics, not backprop, which is why nets still "learned something" and this survived. The book reference (NeuroNet_DNG) fixed this in a later article version; our kernel descended from the earlier one. Confounds every model-based negative verdict to date. Also in this commit, same root cause family: - per-sample UpdateWeightsAdam (OpenCL): input for slot group j was read at matrix_i[j] instead of matrix_i[j*4] (corrupted outer product past group 0), and dispatch dim 1 sized on ceil(inputs/4) left the bias column unreachable whenever inputs%4==0 - dense biases never trained on OpenCL. Rewritten as a lane-guarded scalar loop keeping our Adam conventions (sqrt-stored v, decoupled decay, both clamps, no sign gate). The batched accum path never had either bug; this kernel is what SetBatchSize(1) runs - including online continual learning on client machines, where OpenCL is the only tier. - conv backward passed raw (int)Activation() where the kernels expect NativeActivationCode(): NONE took the tanh branch (clamping a BN layer's unbounded z-scores), TANH took sigmoid, PRELU took none. Dormant only because the conv sits at layer 1 today. - hidden-gradient dispatch over Neurons()+1 dropped to Neurons(): biases get no backprop gradient and the extra work-item only ever read past matrix_o. All three backends (Network.cl, WarriorCPU.cpp, WarriorDML.cpp HLSL) changed in lockstep; DML gained an `inputs` constant to derive the row stride. New dense_backprop_check.cpp proves the CPU kernel is central-finite-difference consistent with the real forward kernel on 8x8, 64x3, 33x64, 5x3 (max diff 3e-9) and that all three activation branches match transcription. All 16 checks pass. Offline math check only - the in-situ proof remains the per-layer dW/W report on a real era. FORCES FULL RETRAIN. Both DLLs rebuilt and redeployed to MQL5\Libraries. Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
1134 lines
60 KiB
Common Lisp
1134 lines
60 KiB
Common Lisp
// fp32 throughout - this project targets consumer GPUs (e.g. AMD Polaris/RX 580) where double
|
|
// precision runs at a small fraction of float throughput. No kernel below needs fp64 range or
|
|
// precision: weight/delta clamps stay within +-100, gradients/activations are float32-safe, and
|
|
// the Adam epsilon guard (vt > 0 ? vt : l*10) isn't a tiny-epsilon comparison sensitive to
|
|
// fp32 vs fp64 underflow. Removing cl_khr_fp64 means any accidental double slipping back in fails
|
|
// kernel compilation loudly instead of silently running at 1/16th speed.
|
|
// Tightened from 1.0e6 - see AI\Network.mqh's matching MAX_WEIGHT comment for the full rationale
|
|
// (that ceiling was so loose it never engaged before training had already gone unstable). 100.0
|
|
// matches the equivalent clamp in Dmitriy Gizlyk's reference NeuroNet.cl engine.
|
|
#define MAX_WEIGHT 100.0f
|
|
// Floor on the magnitude of a saturated tanh/sigmoid unit's activation derivative - see AI\Network.mqh's
|
|
// matching MIN_ACTIVATION_DERIVATIVE comment for the full rationale (a unit pinned near its activation
|
|
// extremes must still receive some corrective gradient, however small, rather than exactly zero).
|
|
// 2026-07-27: Increased from 1.0e-4f to 1.0e-3f to strengthen the escape signal through saturated hidden
|
|
// neurons. On the fp32 OpenCL backend, a saturated sigmoid output layer (the "neutral collapse from era 0"
|
|
// bug) propagated only this floor value backward through every hidden layer, which at 1e-4 was too weak
|
|
// to ever pull weights out of saturation. The output-layer saturation root cause is now fixed (NONE/logit
|
|
// activation instead of SIGMOID — see ExpertSignalAIBase.mqh's BuildFreshTopology), but this increase
|
|
// provides a 10× stronger safety net against any remaining fp32-specific saturation in hidden layers.
|
|
#define MIN_ACTIVATION_DERIVATIVE 1.0e-3f
|
|
// Decoupled (AdamW-style) weight decay applied inside every Adam kernel below - see WarriorCPU.cpp's
|
|
// matching WEIGHT_DECAY comment for the full rationale (unbounded slow weight growth across hundreds
|
|
// of oversampled training eras causing periodic collapse-then-recover cycles). 0.01 is the standard
|
|
// AdamW default.
|
|
#define WEIGHT_DECAY 0.001f
|
|
// Per-step update clip - see WarriorCPU.cpp's matching MAX_WEIGHT_DELTA comment for the full
|
|
// rationale (weight decay alone didn't stop the collapse cycles - they turned out to be sudden Adam
|
|
// overshoot events, most likely from 5x back-to-back oversampling replay building artificially
|
|
// correlated momentum). Applied to the raw delta BEFORE it's added to the weight, unlike MAX_WEIGHT.
|
|
#define MAX_WEIGHT_DELTA 0.1f
|
|
//+------------------------------------------------------------------+
|
|
//| |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void FeedForward(__global float *matrix_w,
|
|
__global float *matrix_i,
|
|
__global float *matrix_o,
|
|
int inputs, int activation)
|
|
{
|
|
int i = get_global_id(0);
|
|
float sum = 0.0f;
|
|
float4 inp, weight;
|
|
int shift = (inputs + 1) * i;
|
|
for(int k = 0; k <= inputs; k = k + 4)
|
|
{
|
|
switch(inputs - k)
|
|
{
|
|
case 0:
|
|
inp = (float4)(1, 0, 0, 0);
|
|
weight = (float4)(matrix_w[shift + k], 0, 0, 0);
|
|
break;
|
|
case 1:
|
|
inp = (float4)(matrix_i[k], 1, 0, 0);
|
|
weight = (float4)(matrix_w[shift + k], matrix_w[shift + k + 1], 0, 0);
|
|
break;
|
|
case 2:
|
|
inp = (float4)(matrix_i[k], matrix_i[k + 1], 1, 0);
|
|
weight = (float4)(matrix_w[shift + k], matrix_w[shift + k + 1], matrix_w[shift + k + 2], 0);
|
|
break;
|
|
case 3:
|
|
inp = (float4)(matrix_i[k], matrix_i[k + 1], matrix_i[k + 2], 1);
|
|
weight = (float4)(matrix_w[shift + k], matrix_w[shift + k + 1], matrix_w[shift + k + 2], matrix_w[shift + k + 3]);
|
|
break;
|
|
default:
|
|
inp = (float4)(matrix_i[k], matrix_i[k + 1], matrix_i[k + 2], matrix_i[k + 3]);
|
|
weight = (float4)(matrix_w[shift + k], matrix_w[shift + k + 1], matrix_w[shift + k + 2], matrix_w[shift + k + 3]);
|
|
break;
|
|
}
|
|
sum += dot(inp, weight);
|
|
}
|
|
switch(activation)
|
|
{
|
|
case 0:
|
|
sum = tanh(sum);
|
|
break;
|
|
case 1:
|
|
sum = 1 / (1 + exp(-clamp(sum, -50.0f, 50.0f)));
|
|
break;
|
|
case 2: // PReLU, param=0.01 - matches FeedForwardConv's case 2 below
|
|
if(sum < 0)
|
|
sum *= 0.01f;
|
|
break;
|
|
}
|
|
matrix_o[i] = sum;
|
|
barrier(CLK_GLOBAL_MEM_FENCE);
|
|
}
|
|
//+------------------------------------------------------------------+
|
|
//| |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void CaclOutputGradient(__global float *matrix_t,
|
|
__global float *matrix_o,
|
|
__global float *matrix_ig,
|
|
int activation)
|
|
{
|
|
int i = get_global_id(0);
|
|
float temp = 0;
|
|
float out = matrix_o[i];
|
|
switch(activation)
|
|
{
|
|
case 0:
|
|
// Deliberately NOT multiplied by the tanh derivative (1-out^2): that
|
|
// factor vanishes as out approaches +-1, which is exactly where a
|
|
// TANH output neuron needs to converge for a +-1 target (e.g. the
|
|
// buy/sell extremes of a single-neuron regression head), stalling
|
|
// training right when it matters most. Using the raw (target-out)
|
|
// error here is the same fix as pairing sigmoid with cross-entropy -
|
|
// a bounded activation whose output error isn't re-damped by its
|
|
// own saturating derivative.
|
|
temp = clamp(matrix_t[i], -1.0f, 1.0f) - out;
|
|
break;
|
|
case 1:
|
|
// Also deliberately NOT multiplied by the sigmoid derivative out*(1-out) - same
|
|
// reasoning as case 0 above. This is the classification output layer (3 neurons,
|
|
// one-hot 0/1 targets - see BuildFreshTopology()); (target-out) unscaled is exactly
|
|
// the sigmoid+binary-cross-entropy gradient (the loss derivative and the sigmoid
|
|
// derivative algebraically cancel), whereas re-multiplying by out*(1-out) here is
|
|
// the MSE-with-sigmoid formula, which damps toward zero as out approaches 0 or 1 -
|
|
// exactly where a 0/1-target neuron needs to converge. With this damping (plus the
|
|
// 0.00000001 underflow guard at the saturated extremes), the classification output
|
|
// layer's effective learning rate was throttled to a small fraction of the
|
|
// regression path's, observed in practice as all 3 output neurons converging to
|
|
// (and never escaping) an identical value - softmax of three equal logits is always
|
|
// exactly 1/3 each, with ties resolved to index 0 ("Buy" wins by default).
|
|
temp = clamp(matrix_t[i], 0.0f, 1.0f) - out;
|
|
break;
|
|
case 2:
|
|
// PReLU, param=0.01 - unlike tanh/sigmoid, PReLU is unbounded and never saturates,
|
|
// so none of the clamp-to-range/derivative-skip workarounds above are needed here:
|
|
// plain (target-out) scaled by the derivative at this output is exact.
|
|
temp = (matrix_t[i] - out) * (out >= 0 ? 1.0f : 0.01f);
|
|
break;
|
|
default:
|
|
// NONE (raw logits - see the 3-class softmax classification output layer in
|
|
// BuildFreshTopology()/NativeActivationCode()'s NONE case, which maps here): identity
|
|
// activation, derivative 1 everywhere, so the plain (target-out) error passes through
|
|
// unscaled - same formula CNeuron::calcOutputGradients() already uses unconditionally
|
|
// on the un-accelerated CPU fallback path. Leaving this case unhandled left `temp` at
|
|
// its 0.0 initializer, i.e. every softmax-classification output neuron got a
|
|
// permanently zero gradient on this backend - which, since calcHiddenGradients()
|
|
// propagates backward FROM the output layer's gradient, silently froze the entire
|
|
// network (not just the output layer) at its random initial weights whenever an
|
|
// OpenCL-accelerated tier was active.
|
|
temp = matrix_t[i] - out;
|
|
break;
|
|
}
|
|
matrix_ig[i] = temp;
|
|
}
|
|
//+------------------------------------------------------------------+
|
|
//| |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void CaclHiddenGradient(__global float *matrix_w,
|
|
__global float *matrix_g,
|
|
__global float *matrix_o,
|
|
__global float *matrix_ig,
|
|
int outputs, int activation)
|
|
{
|
|
int i = get_global_id(0);
|
|
// This layer's neuron count. Dispatch is exactly Neurons() (no +1: biases receive no backprop
|
|
// gradient, and the old +1 work-item only ever produced an out-of-bounds read of matrix_o).
|
|
int inputs = get_global_size(0);
|
|
float sum = 0;
|
|
float out = matrix_o[i];
|
|
// matrix_w is THIS layer's outgoing weight matrix, laid out exactly as FeedForward consumes it:
|
|
// one row per NEXT-layer neuron k, row stride (inputs + 1), so the weight FROM this neuron i
|
|
// INTO next-layer neuron k is matrix_w[k * (inputs + 1) + i] - a column read, dL/dout_i =
|
|
// sum_k g[k] * W[k][i]. Until 2026-08-11 this kernel read matrix_w[(outputs + 1) * i + k]: the
|
|
// TRANSPOSE for square layers and, for the non-square boundaries this EA actually builds
|
|
// (tapered stacks, the 3-neuron head), a mis-strided walk that ran past the buffer's end. Every
|
|
// gradient crossing a dense boundary on its way down - the entire learning signal reaching the
|
|
// BN/conv/LSTM front ends - passed through that fixed wrong matrix, i.e. feedback-alignment
|
|
// dynamics rather than backprop (which is why nets still "learned something" and it survived).
|
|
// The book reference (references\...\NeuroNet_DNG\NeuroNet.cl CalcHiddenGradient) has the
|
|
// corrected form; this repo's kernel descended from the earlier article version. Mirrored in
|
|
// DirectML\WarriorCPU.cpp (CPU_CalcHiddenGradient) and WarriorDML.cpp (kHlslHiddenGradient) -
|
|
// the three must stay in lockstep.
|
|
for(int k = 0; k < outputs; k++)
|
|
sum += matrix_g[k] * matrix_w[k * (inputs + 1) + i];
|
|
switch(activation)
|
|
{
|
|
case 0:
|
|
sum = clamp(sum + out, -1.0f, 1.0f) - out;
|
|
sum = sum * fmax(MIN_ACTIVATION_DERIVATIVE, 1 - pow(out, 2.0f));
|
|
break;
|
|
case 1:
|
|
sum = clamp(sum + out, 0.0f, 1.0f) - out;
|
|
sum = sum * fmax(MIN_ACTIVATION_DERIVATIVE, out * (1 - out));
|
|
break;
|
|
case 2:
|
|
// PReLU, param=0.01 - unbounded, non-saturating, so the plain chain rule applies
|
|
// directly with no clamp/implied-target reformulation needed (contrast with the
|
|
// tanh/sigmoid cases above, which exist specifically to dodge saturation).
|
|
sum = sum * (out >= 0 ? 1.0f : 0.01f);
|
|
break;
|
|
}
|
|
matrix_ig[i] = sum;
|
|
}
|
|
//+------------------------------------------------------------------+
|
|
//| |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void UpdateWeightsMomentum(__global float *matrix_w,
|
|
__global float *matrix_g,
|
|
__global float *matrix_i,
|
|
__global float *matrix_dw,
|
|
int inputs, float learning_rates, float momentum, int optimizer)
|
|
{
|
|
int i = get_global_id(0);
|
|
int j = get_global_id(1);
|
|
int wi = i * (inputs + 1) + j;
|
|
// `optimizer` is retained (unused) only to keep this kernel's argument list binary-compatible with
|
|
// the already-deployed WarriorDML.dll/WarriorCPU.dll exports. It used to select an index-parity sign
|
|
// flip on the gradient ("DFA"), which put half of every weight tensor into permanent gradient ASCENT
|
|
// - see ENUM_OPTIMIZATION's comment in AI\Network.mqh. Always plain descent now; callers pass 0.
|
|
float delta = learning_rates * (matrix_g[i] * (j < inputs ? matrix_i[j] : 1)) + momentum * matrix_dw[wi];
|
|
matrix_dw[wi] = delta;
|
|
// Unlike UpdateWeightsConvMomentum/UpdateWeightsConvAdam/LSTM_UpdateWeightsAdam below, this dense
|
|
// kernel (and UpdateWeightsAdam right after it) had never gotten the MAX_WEIGHT clamp - a gradient
|
|
// spike could drive a weight to +-Infinity with nothing to stop it, and the next Adam step turns
|
|
// Infinity into NaN (Inf/Inf), which then propagates through every FeedForward sum touching that
|
|
// weight and never recovers. Dense layers are exactly where the classification output head lives,
|
|
// so this is the most consequential of the three missing clamps, not the least.
|
|
matrix_w[wi] = clamp(matrix_w[wi] + delta, -MAX_WEIGHT, MAX_WEIGHT);
|
|
};
|
|
//+------------------------------------------------------------------+
|
|
//| |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void UpdateWeightsAdam(__global float *matrix_w,
|
|
__global const float *matrix_g,
|
|
__global const float *matrix_i,
|
|
__global float *matrix_m,
|
|
__global float *matrix_v,
|
|
const int inputs, const float l, const float b1, const float b2)
|
|
{
|
|
const int i = get_global_id(0);
|
|
const int j = get_global_id(1);
|
|
// Each work-item owns up to 4 consecutive slots of weight row i. The row is (inputs + 1) wide -
|
|
// the last column is the bias, whose implicit input is 1 - and dispatch dim 1 is sized on
|
|
// ceil((inputs + 1) / 4), so the bias column is always reachable. Until 2026-08-11 this kernel
|
|
// had two defects the book reference had already fixed: the input for slot group j was read at
|
|
// matrix_i[j] instead of matrix_i[j * 4] (every group past the first paired its weights with the
|
|
// WRONG input - a corrupted outer product), and the tail switch tested (inputs - j * 4) with
|
|
// dim 1 sized on ceil(inputs / 4), which made the bias column unreachable whenever
|
|
// inputs % 4 == 0 (32/64-wide layers: dense biases simply never trained on OpenCL). The batched
|
|
// path (AccumulateWeightGrad + ApplyAccumAdam) never had either bug; this per-sample kernel is
|
|
// what SetBatchSize(1) paths run - including online continual learning on client machines,
|
|
// where OpenCL is the only tier. The per-lane guard keeps the tail group off the next row.
|
|
for(int n = 0; n < 4; n++)
|
|
{
|
|
int col = j * 4 + n;
|
|
if(col > inputs)
|
|
break;
|
|
int wi = i * (inputs + 1) + col;
|
|
float inp = (col < inputs ? matrix_i[col] : 1.0f);
|
|
float g = matrix_g[i] * inp;
|
|
float mt = b1 * matrix_m[wi] + (1 - b1) * g;
|
|
// v is STORED already square-rooted (a standard deviation) so it can be the denominator
|
|
// directly; it must therefore be squared back before re-entering the recursion. Feeding the
|
|
// stored sqrt in as if it were the variance - which every Adam kernel here did until
|
|
// 2026-08-09 - mixes a std dev with a variance, and the resulting recursion has a fixed
|
|
// point at v ~= b2 for ANY |g| below 1. The denominator then stops tracking the gradient
|
|
// scale, delta becomes proportional to |g|, and Adam silently degrades to plain SGD. That is
|
|
// not academic: it is why the conv/LSTM front end (which sits behind a batch-norm and so
|
|
// receives gradients divided by sqrt(var) ~ 500) froze at dW/W 0.00% while the dense stack
|
|
// trained. Measured with DirectML\batch_accum_check.cpp against the shipped kernel: 3285x
|
|
// less displacement at |g|=1e-5 than at |g|=1, where a scale-invariant optimizer gives the
|
|
// same distance for both. See NeuronBatchNorm.mqh, which got this right first.
|
|
float vt = sqrt(b2 * matrix_v[wi] * matrix_v[wi] + (1 - b2) * g * g);
|
|
float w = matrix_w[wi];
|
|
float delta = clamp(l * mt / (vt > 0 ? vt : l * 10) - l * WEIGHT_DECAY * w, -MAX_WEIGHT_DELTA, MAX_WEIGHT_DELTA);
|
|
// See UpdateWeightsMomentum's comment just above for why the matrix_w write is clamped to
|
|
// MAX_WEIGHT - this dense Adam kernel is the one that updates the classification output
|
|
// layer. No sign-agreement gate on the matrix_w write (removed 2026-07, all four backends in
|
|
// sync) - see AI\NeuronCPU.mqh's CNeuron::updateInputWeights comment for why: it rectified
|
|
// the one-hot softmax-CCE gradient stream into a permanent downward ratchet (the all-Neutral
|
|
// collapse).
|
|
matrix_w[wi] = clamp(w + delta, -MAX_WEIGHT, MAX_WEIGHT);
|
|
matrix_m[wi] = mt;
|
|
matrix_v[wi] = vt;
|
|
}
|
|
};
|
|
//+------------------------------------------------------------------+
|
|
//| Max-pooling layer (CNeuronPoolOCL) - no weights, just a sliding |
|
|
//| max over the previous layer's output. Ported from the NeuroNet_DNG|
|
|
//| reference's CNeuronProofOCL kernels (these two had no bug). |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void FeedForwardProof(__global float *matrix_i,
|
|
__global float *matrix_o,
|
|
int inputs, int window, int step)
|
|
{
|
|
int i = get_global_id(0);
|
|
int pos = i * step;
|
|
// Mirrors the bounds guard in this kernel's pure-MQL5 host counterpart
|
|
// (CNeuronPoolOCL::feedForwardCPU in AI\NeuronOCLConvPool.mqh) - outputs is sized so pos shouldn't
|
|
// reach inputs in practice, but without this guard an out-of-bounds pos would read past matrix_i.
|
|
if(pos >= inputs)
|
|
return;
|
|
float result = matrix_i[pos];
|
|
for(int k = 1; k < window; k++)
|
|
{
|
|
int shift = k + pos;
|
|
if(shift >= inputs)
|
|
break;
|
|
result = max(result, matrix_i[shift]);
|
|
}
|
|
matrix_o[i] = result;
|
|
}
|
|
//+------------------------------------------------------------------+
|
|
//| |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void CalcInputGradientProof(__global float *matrix_i,
|
|
__global float *matrix_g,
|
|
__global float *matrix_o,
|
|
__global float *matrix_ig,
|
|
int outputs, int window, int step)
|
|
{
|
|
int i = get_global_id(0);
|
|
float prev_gradient = 0.0f;
|
|
float value = matrix_i[i];
|
|
int start = i - window + step;
|
|
start = (start - start % step) / step;
|
|
int stop = (i - i % step) / step + 1;
|
|
for(int out = max(0, start); out < min(outputs, stop); out++)
|
|
if(value == matrix_o[out])
|
|
prev_gradient += matrix_g[out];
|
|
matrix_ig[i] = prev_gradient;
|
|
}
|
|
//+------------------------------------------------------------------+
|
|
//| Convolution layer (CNeuronConvOCL) - ported from the NeuroNet_DNG |
|
|
//| reference library and adapted to this project's buffer convention.|
|
|
//| window_in+1 weights (incl. bias) are shared across all sliding |
|
|
//| positions; window_out lets one position produce several output |
|
|
//| channels (filters) at once. |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void FeedForwardConv(__global float *matrix_w,
|
|
__global float *matrix_i,
|
|
__global float *matrix_o,
|
|
int inputs, int step, int window_in, int window_out,
|
|
int activation)
|
|
{
|
|
int i = get_global_id(0);
|
|
int shift_out = window_out * i;
|
|
int shift_in = step * i;
|
|
for(int out = 0; out < window_out; out++)
|
|
{
|
|
int shift = (window_in + 1) * out;
|
|
int stop = (window_in <= (inputs - shift_in) ? window_in : (inputs - shift_in));
|
|
float sum = 0.0f;
|
|
for(int k = 0; k < stop; k++)
|
|
sum += matrix_i[shift_in + k] * matrix_w[shift + k];
|
|
sum += matrix_w[shift + window_in]; // bias
|
|
switch(activation)
|
|
{
|
|
case 0:
|
|
sum = tanh(sum);
|
|
break;
|
|
case 1:
|
|
sum = 1 / (1 + exp(-clamp(sum, -50.0f, 50.0f)));
|
|
break;
|
|
case 2: // PReLU, param=0.01
|
|
if(sum < 0)
|
|
sum *= 0.01f;
|
|
break;
|
|
}
|
|
matrix_o[out + shift_out] = sum;
|
|
}
|
|
}
|
|
//+------------------------------------------------------------------+
|
|
//| |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void CalcHiddenGradientConv(__global float *matrix_w,
|
|
__global float *matrix_g,
|
|
__global float *matrix_o,
|
|
__global float *matrix_ig,
|
|
int outputs, int step, int window_in, int window_out,
|
|
int activation)
|
|
{
|
|
int i = get_global_id(0);
|
|
float out = matrix_o[i];
|
|
int start = i - window_in + step;
|
|
start = max((start - start % step) / step, 0);
|
|
int stop = (i - i % step) / step + 1;
|
|
if(stop > (outputs / window_out))
|
|
stop = outputs / window_out;
|
|
float sum = 0.0f;
|
|
for(int h = 0; h < window_out; h++)
|
|
{
|
|
for(int k = start; k < stop; k++)
|
|
{
|
|
int shift_w = (stop - k - 1) * step + i % step + h * (window_in + 1);
|
|
int shift_g = k * window_out + h;
|
|
if(shift_g >= outputs || shift_w >= (window_in + 1) * window_out)
|
|
break;
|
|
sum += matrix_g[shift_g] * matrix_w[shift_w];
|
|
}
|
|
}
|
|
switch(activation)
|
|
{
|
|
case 0:
|
|
sum = clamp(sum + out, -1.0f, 1.0f) - out;
|
|
sum = sum * fmax(MIN_ACTIVATION_DERIVATIVE, 1 - pow(out, 2.0f));
|
|
break;
|
|
case 1:
|
|
sum = clamp(sum + out, 0.0f, 1.0f) - out;
|
|
sum = sum * fmax(MIN_ACTIVATION_DERIVATIVE, out * (1 - out));
|
|
break;
|
|
case 2: // PReLU, param=0.01
|
|
if(out < 0)
|
|
sum *= 0.01f;
|
|
break;
|
|
}
|
|
matrix_ig[i] = sum;
|
|
}
|
|
//+------------------------------------------------------------------+
|
|
//| |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void UpdateWeightsConvMomentum(__global float *matrix_w,
|
|
__global float *matrix_g,
|
|
__global float *matrix_i,
|
|
__global float *matrix_dw,
|
|
int inputs, float learning_rates, float momentum,
|
|
int window_in, int window_out, int step, int optimizer)
|
|
{
|
|
int i = get_global_id(0);
|
|
int shift = i % (window_in + 1);
|
|
int shift_out = (i - shift) / (window_in + 1);
|
|
int total = (inputs - window_in) % step;
|
|
total = (inputs - window_in - total) / step + (total > 0 ? 1 : 0);
|
|
float grad = 0.0f;
|
|
for(int t = 0; t < total; t++)
|
|
{
|
|
if(shift != window_in && (shift + t * step) >= inputs)
|
|
break;
|
|
grad += matrix_g[t * window_out + shift_out] * (shift == window_in ? 1 : matrix_i[shift + t * step]);
|
|
}
|
|
// `optimizer` retained unused - see UpdateWeightsMomentum's comment above.
|
|
float delta = learning_rates * grad + momentum * matrix_dw[i];
|
|
matrix_dw[i] = delta;
|
|
matrix_w[i] = clamp(matrix_w[i] + delta, -MAX_WEIGHT, MAX_WEIGHT);
|
|
};
|
|
//+------------------------------------------------------------------+
|
|
//| |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void UpdateWeightsConvAdam(__global float *matrix_w,
|
|
__global const float *matrix_g,
|
|
__global const float *matrix_i,
|
|
__global float *matrix_m,
|
|
__global float *matrix_v,
|
|
const int inputs, const float l, const float b1, const float b2,
|
|
int window_in, int window_out, int step)
|
|
{
|
|
int i = get_global_id(0); // weight offset within a single output channel's (window_in+1) block
|
|
if(i > window_in)
|
|
return;
|
|
int total = (inputs - (window_in - step)) % step;
|
|
total = (inputs - (window_in - step) - total) / step + (total > 0 ? 1 : 0);
|
|
for(int out = 0; out < window_out; out++)
|
|
{
|
|
int shift_w = i + out * (window_in + 1);
|
|
float grad = 0.0f;
|
|
for(int t = 0; t < total; t++)
|
|
{
|
|
if(i != window_in && (i + t * step) >= inputs)
|
|
break;
|
|
grad += matrix_g[t * window_out + out] * (i == window_in ? 1 : matrix_i[i + t * step]);
|
|
}
|
|
float mt = b1 * matrix_m[shift_w] + (1 - b1) * grad;
|
|
// v squared back before the recursion - see UpdateWeightsAdam's comment above for why.
|
|
float vt = sqrt(b2 * matrix_v[shift_w] * matrix_v[shift_w] + (1 - b2) * pow(grad, 2.0f));
|
|
float delta = clamp(l * mt / (vt > 0 ? vt : l * 10) - l * WEIGHT_DECAY * matrix_w[shift_w], -MAX_WEIGHT_DELTA, MAX_WEIGHT_DELTA);
|
|
// Sign-agreement gate removed - see UpdateWeightsAdam's comment above for why.
|
|
matrix_w[shift_w] = clamp(matrix_w[shift_w] + delta, -MAX_WEIGHT, MAX_WEIGHT);
|
|
matrix_m[shift_w] = mt;
|
|
matrix_v[shift_w] = vt;
|
|
}
|
|
};
|
|
//+------------------------------------------------------------------+
|
|
//| MINI-BATCH GRADIENT ACCUMULATION (2026-08-09 audit, F4). |
|
|
//| |
|
|
//| Training here is otherwise pure online SGD - one weight update |
|
|
//| per sample - which makes the per-update gradient noise maximal |
|
|
//| and the end-of-era weight state a high-variance random draw. The |
|
|
//| fix is to sum each weight's gradient over TRAIN_BATCH_SIZE |
|
|
//| samples and take ONE optimizer step on the mean. |
|
|
//| |
|
|
//| Only the ACCUMULATION lives here. The optimizer step itself is |
|
|
//| elementwise O(weights), so - exactly like the batch-norm layer |
|
|
//| and the softmax+CCE output gradient before it - it is computed |
|
|
//| once per BATCH in MQL5 against the host mirrors of these buffers |
|
|
//| (CNeuronBaseOCL::ApplyAccumulatedGradients). That keeps ONE |
|
|
//| optimizer implementation for all four compute tiers instead of |
|
|
//| four that can drift apart, which this engine has been bitten by |
|
|
//| repeatedly. |
|
|
//| |
|
|
//| The per-weight gradient expressions below are copied verbatim |
|
|
//| from UpdateWeightsAdam / UpdateWeightsConvAdam above - they MUST |
|
|
//| stay identical, because batch size 1 has to reproduce the |
|
|
//| unbatched path exactly. |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void AccumulateWeightGrad(__global float *matrix_acc,
|
|
__global const float *matrix_g,
|
|
__global const float *matrix_i,
|
|
const int inputs)
|
|
{
|
|
const int i = get_global_id(0); // destination neuron
|
|
const int j = get_global_id(1); // source slot, j == inputs is the bias
|
|
const int wi = i * (inputs + 1) + j;
|
|
// Deliberately scalar, not the float4 form UpdateWeightsAdam uses: this kernel replaces that
|
|
// dispatch on batched runs rather than adding to it, so the dispatch count is unchanged, and a
|
|
// plain expression is far easier to hold against the host-side apply it has to agree with.
|
|
matrix_acc[wi] += matrix_g[i] * (j < inputs ? matrix_i[j] : 1.0f);
|
|
};
|
|
//+------------------------------------------------------------------+
|
|
//| Convolution counterpart. Flat over (window_in+1)*window_out, i.e. |
|
|
//| UpdateWeightsConvMomentum's indexing, but with the POSITION COUNT |
|
|
//| taken from UpdateWeightsConvAdam. |
|
|
//| |
|
|
//| Those two disagree, and it is worth being explicit about which is |
|
|
//| reproduced here. Positions of a window_in-wide window stepping by |
|
|
//| `step` over `inputs` values is floor((inputs-window_in)/step)+1. |
|
|
//| The Adam kernel's total - ceil((inputs-window_in+step)/step) - |
|
|
//| equals that; the Momentum kernel's ceil((inputs-window_in)/step) |
|
|
//| is one SHORT whenever step divides the span, so the SGD conv path |
|
|
//| silently drops the last sliding position's gradient. That is a |
|
|
//| pre-existing defect in the SGD path (ADAM is the default and the |
|
|
//| only optimizer these models ship with) and is NOT fixed here - |
|
|
//| fixing it would change SGD's unbatched behaviour, which is not |
|
|
//| this change's business. The correct count is used because this |
|
|
//| accumulator feeds the shared host-side apply for BOTH optimizers. |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void AccumulateWeightGradConv(__global float *matrix_acc,
|
|
__global const float *matrix_g,
|
|
__global const float *matrix_i,
|
|
const int inputs, int window_in, int window_out, int step)
|
|
{
|
|
int i = get_global_id(0);
|
|
int shift = i % (window_in + 1);
|
|
int shift_out = (i - shift) / (window_in + 1);
|
|
int total = (inputs - (window_in - step)) % step;
|
|
total = (inputs - (window_in - step) - total) / step + (total > 0 ? 1 : 0);
|
|
float grad = 0.0f;
|
|
for(int t = 0; t < total; t++)
|
|
{
|
|
if(shift != window_in && (shift + t * step) >= inputs)
|
|
break;
|
|
grad += matrix_g[t * window_out + shift_out] * (shift == window_in ? 1 : matrix_i[shift + t * step]);
|
|
}
|
|
matrix_acc[i] += grad;
|
|
};
|
|
//+------------------------------------------------------------------+
|
|
//| dst += src, elementwise. Used by the LSTM's mini-batch path: that |
|
|
//| layer's WeightsGradient already holds one sample's complete dW |
|
|
//| (summed over the BPTT timesteps), but the seq-backward kernels |
|
|
//| MEMSET it at the top of every call, so the running batch total |
|
|
//| has to live in a second buffer. Deliberately generic rather than |
|
|
//| LSTM-specific - it is just a vector add. |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void AccumulateBufferInto(__global float *dst,
|
|
__global const float *src)
|
|
{
|
|
const int i = get_global_id(0);
|
|
dst[i] += src[i];
|
|
};
|
|
//+------------------------------------------------------------------+
|
|
//| MINI-BATCH APPLY - one optimizer step on the accumulated mean. |
|
|
//| |
|
|
//| These replace a HOST-SIDE apply that read the whole weight matrix,|
|
|
//| the accumulator and both Adam moments back over the bus, stepped |
|
|
//| them in MQL5, and wrote four buffers out again: eight full-matrix |
|
|
//| transfers per batch per weight block, each a blocking sync. At |
|
|
//| batch 8 that was ~1 weight-matrix-worth of PCIe traffic PER |
|
|
//| SAMPLE, which is why an RX 580 lost to a CPU thread pool on this |
|
|
//| workload. Market builds cannot import a DLL, so OpenCL is the |
|
|
//| tier clients actually run and this path is the product's speed. |
|
|
//| |
|
|
//| Flat and elementwise: the accumulator already holds the summed |
|
|
//| per-weight gradient (AccumulateWeightGrad above), so there is no |
|
|
//| outer product to reconstruct and no float4 dispatch to mirror. |
|
|
//| |
|
|
//| `scale` is 1/batchCount - the MEAN. The arithmetic is a |
|
|
//| line-for-line transcription of UpdateWeightsAdam above, INCLUDING |
|
|
//| the squared-back second moment, both clamps and the decoupled |
|
|
//| decay, so batch size 1 still reproduces the unbatched kernel. The |
|
|
//| host implementation in NeuronOCLBase.mqh remains the reference |
|
|
//| for the DLL and pure-MQL5 tiers and must be edited in step. |
|
|
//| |
|
|
//| The accumulator is zeroed HERE, saving a separate clear dispatch |
|
|
//| and guaranteeing it cannot be left dirty for the next batch by an |
|
|
//| early return between the step and the clear. |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void ApplyAccumAdam(__global float *matrix_w,
|
|
__global float *matrix_acc,
|
|
__global float *matrix_m,
|
|
__global float *matrix_v,
|
|
const float scale, const float l,
|
|
const float b1, const float b2)
|
|
{
|
|
const int i = get_global_id(0);
|
|
const float grad = matrix_acc[i] * scale;
|
|
const float wv = matrix_w[i];
|
|
const float mt = b1 * matrix_m[i] + (1 - b1) * grad;
|
|
// Squared back before re-entering the recursion - see UpdateWeightsAdam for the full account of
|
|
// what feeding the stored square root in as the variance did to this engine.
|
|
const float vt = sqrt(b2 * matrix_v[i] * matrix_v[i] + (1 - b2) * grad * grad);
|
|
const float delta = clamp(l * mt / (vt > 0 ? vt : l * 10) - l * WEIGHT_DECAY * wv,
|
|
-MAX_WEIGHT_DELTA, MAX_WEIGHT_DELTA);
|
|
matrix_w[i] = clamp(wv + delta, -MAX_WEIGHT, MAX_WEIGHT);
|
|
matrix_m[i] = mt;
|
|
matrix_v[i] = vt;
|
|
matrix_acc[i] = 0.0f;
|
|
};
|
|
//+------------------------------------------------------------------+
|
|
//| SGD/momentum counterpart. Note the delta is NOT clamped to |
|
|
//| MAX_WEIGHT_DELTA here - matching UpdateWeightsMomentum and the |
|
|
//| host apply, both of which clamp only the resulting weight. Kept |
|
|
//| that way deliberately rather than "improved": the batched path |
|
|
//| has to reproduce the unbatched one, not diverge from it. |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void ApplyAccumMomentum(__global float *matrix_w,
|
|
__global float *matrix_acc,
|
|
__global float *matrix_dw,
|
|
const float scale, const float lr, const float momentum)
|
|
{
|
|
const int i = get_global_id(0);
|
|
const float grad = matrix_acc[i] * scale;
|
|
const float delta = lr * grad + momentum * matrix_dw[i];
|
|
matrix_dw[i] = delta;
|
|
matrix_w[i] = clamp(matrix_w[i] + delta, -MAX_WEIGHT, MAX_WEIGHT);
|
|
matrix_acc[i] = 0.0f;
|
|
};
|
|
//+------------------------------------------------------------------+
|
|
//| BATCH NORM (CNeuronBatchNormOCL), device-side since 2026-08-09. |
|
|
//| |
|
|
//| This layer computed host-side on every tier, which on OpenCL |
|
|
//| meant a blocking read of the previous layer's output and a |
|
|
//| blocking write of its own - per sample, per layer, forward AND |
|
|
//| backward. These kernels are LINE-FOR-LINE transcriptions of the |
|
|
//| host implementations in AI\NeuronBatchNorm.mqh (NormalizeHost, |
|
|
//| calcInputGradients, StepGammaBeta), including every NaN guard and |
|
|
//| clamp - that layer latched NaN into its persisted statistics once |
|
|
//| (2026-08-02) and the guards are what keep a single garbage value |
|
|
//| from poisoning the .nnw forever. The host copies remain the |
|
|
//| reference AND the runtime for the DLL / pure-MQL5 tiers; edit the |
|
|
//| two together. The MQL5 side also self-checks each kernel against |
|
|
//| its host twin on first use and latches back to the host path on |
|
|
//| any disagreement, so a transcription error costs a warning and |
|
|
//| some speed, never a corrupted model. |
|
|
//| |
|
|
//| Slot layout below MUST match NeuronBatchNorm.mqh's BN_OPT_*. |
|
|
//+------------------------------------------------------------------+
|
|
#define BN_OPT_STRIDE 9
|
|
#define BN_OPT_MEAN 0
|
|
#define BN_OPT_VAR 1
|
|
#define BN_OPT_NX 2
|
|
#define BN_OPT_GAMMA 3
|
|
#define BN_OPT_BETA 4
|
|
#define BN_OPT_MG 5
|
|
#define BN_OPT_MB 6
|
|
#define BN_OPT_VG 7
|
|
#define BN_OPT_VB 8
|
|
#define BN_EPSILON 1.0e-10f
|
|
#define BN_MIN_STD 1.0e-4f
|
|
#define BN_MAX_INPUT 1.0e6f
|
|
//+------------------------------------------------------------------+
|
|
//| Forward transform for one unit. `w` is the host-computed |
|
|
//| effective EMA window (the 1,2,3..iBatchSize ramp lives host-side |
|
|
//| with iSamplesSeen); `frozen` mirrors bStatsFrozen - statistics |
|
|
//| are USED but not UPDATED, classic inference semantics. |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void BatchNormForward(__global const float *matrix_i,
|
|
__global float *matrix_o,
|
|
__global float *options,
|
|
const float w, const int frozen)
|
|
{
|
|
const int i = get_global_id(0);
|
|
const int shift = i * BN_OPT_STRIDE;
|
|
float x = matrix_i[i];
|
|
// Guards in the same order as NormalizeHost: non-finite -> 0, then the range clamp that keeps the
|
|
// statistics arithmetic bounded (a merely-huge finite x overflows (x-mean)^2 long before infinity).
|
|
if(isnan(x) || isinf(x))
|
|
x = 0.0f;
|
|
x = clamp(x, -BN_MAX_INPUT, BN_MAX_INPUT);
|
|
float mean = options[shift + BN_OPT_MEAN];
|
|
float variance = options[shift + BN_OPT_VAR];
|
|
// Self-heal statistics poisoned before the guards existed - Load faithfully restores NaN.
|
|
if(isnan(mean) || isinf(mean))
|
|
mean = x;
|
|
if(isnan(variance) || isinf(variance) || variance < 0.0f)
|
|
variance = 0.0f;
|
|
if(frozen == 0)
|
|
{
|
|
// EMA update, new mean used by the variance term - exact order of the host code.
|
|
mean = (mean * (w - 1.0f) + x) / w;
|
|
variance = (variance * (w - 1.0f) + (x - mean) * (x - mean)) / w;
|
|
}
|
|
const float delta = x - mean;
|
|
const float sd = fmax(sqrt(variance + BN_EPSILON), BN_MIN_STD);
|
|
const float nx = delta / sd;
|
|
float gamma = options[shift + BN_OPT_GAMMA];
|
|
float beta = options[shift + BN_OPT_BETA];
|
|
// Local substitution only, exactly like the host forward: the PERSISTED repair happens in the
|
|
// gamma/beta step, which is the copy that makes it stick.
|
|
if(isnan(gamma) || isinf(gamma))
|
|
gamma = 1.0f;
|
|
if(isnan(beta) || isinf(beta))
|
|
beta = 0.0f;
|
|
if(frozen == 0)
|
|
{
|
|
options[shift + BN_OPT_MEAN] = mean;
|
|
options[shift + BN_OPT_VAR] = variance;
|
|
}
|
|
// nx cached even when frozen - costs nothing, keeps the buffer consistent with the output.
|
|
options[shift + BN_OPT_NX] = nx;
|
|
matrix_o[i] = gamma * nx + beta;
|
|
};
|
|
//+------------------------------------------------------------------+
|
|
//| Backward for one unit: dL/dx = gamma/sd * dL/dy, then the |
|
|
//| PREVIOUS layer's activation derivative applied with the same |
|
|
//| clamp-to-range "implied target" treatment CaclHiddenGradient |
|
|
//| uses, so from below a batch-norm layer is indistinguishable from |
|
|
//| any other. `activation` is the NATIVE code (0=TANH 1=SIGMOID |
|
|
//| 2=PRELU, anything else passes through) - NativeActivationCode() |
|
|
//| at the dispatch site, never a raw enum cast. |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void BatchNormHiddenGrad(__global const float *matrix_g,
|
|
__global const float *prev_o,
|
|
__global float *prev_g,
|
|
__global const float *options,
|
|
const int activation)
|
|
{
|
|
const int i = get_global_id(0);
|
|
const int shift = i * BN_OPT_STRIDE;
|
|
const float sd = fmax(sqrt(options[shift + BN_OPT_VAR] + BN_EPSILON), BN_MIN_STD);
|
|
float g = matrix_g[i] * options[shift + BN_OPT_GAMMA] / sd;
|
|
const float out = prev_o[i];
|
|
switch(activation)
|
|
{
|
|
case 0:
|
|
g = clamp(g + out, -1.0f, 1.0f) - out;
|
|
g = g * fmax(MIN_ACTIVATION_DERIVATIVE, 1.0f - out * out);
|
|
break;
|
|
case 1:
|
|
g = clamp(g + out, 0.0f, 1.0f) - out;
|
|
g = g * fmax(MIN_ACTIVATION_DERIVATIVE, out * (1.0f - out));
|
|
break;
|
|
case 2:
|
|
g = g * (out >= 0.0f ? 1.0f : 0.01f);
|
|
break;
|
|
default:
|
|
break;
|
|
}
|
|
prev_g[i] = g;
|
|
};
|
|
//+------------------------------------------------------------------+
|
|
//| Mini-batch accumulate for gamma/beta: dL/dgamma = dL/dy * nx, |
|
|
//| dL/dbeta = dL/dy, summed per unit into acc[2i]/acc[2i+1]. Formed |
|
|
//| per sample because BN_OPT_NX is overwritten by every forward. |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void BatchNormAccumGammaBeta(__global const float *matrix_g,
|
|
__global const float *options,
|
|
__global float *acc)
|
|
{
|
|
const int i = get_global_id(0);
|
|
const float g = matrix_g[i];
|
|
acc[2 * i] += g * options[i * BN_OPT_STRIDE + BN_OPT_NX];
|
|
acc[2 * i + 1] += g;
|
|
};
|
|
//+------------------------------------------------------------------+
|
|
//| One gamma/beta optimizer step on the accumulated mean, then the |
|
|
//| accumulator is zeroed. Transcribed from StepGammaBeta including |
|
|
//| its exact write ordering: moments are stored BEFORE the |
|
|
//| finiteness check, so a skipped unit still advances its optimizer |
|
|
//| state, and the persisted gamma/beta self-heal happens HERE (this |
|
|
//| is the copy that sticks). NO WEIGHT DECAY - decaying gamma toward |
|
|
//| zero is the collapse this layer exists to prevent. The per-sample |
|
|
//| path is this same pair of kernels at scale=1. |
|
|
//| `optimizer`: 0 = SGD/momentum (lr=eta, momentum=alpha), else Adam |
|
|
//| (lt = bias-corrected step, b1, b2; v stored square-rooted, same |
|
|
//| convention as every Adam kernel in this file). |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void BatchNormApplyGammaBeta(__global float *options,
|
|
__global float *acc,
|
|
const float scale, const float lt,
|
|
const float b1, const float b2,
|
|
const float lr, const float momentum,
|
|
const int optimizer)
|
|
{
|
|
const int i = get_global_id(0);
|
|
const int shift = i * BN_OPT_STRIDE;
|
|
const float gGamma = acc[2 * i] * scale;
|
|
const float gBeta = acc[2 * i + 1] * scale;
|
|
acc[2 * i] = 0.0f;
|
|
acc[2 * i + 1] = 0.0f;
|
|
float gamma = options[shift + BN_OPT_GAMMA];
|
|
float beta = options[shift + BN_OPT_BETA];
|
|
if(isnan(gamma) || isinf(gamma))
|
|
gamma = 1.0f;
|
|
if(isnan(beta) || isinf(beta))
|
|
beta = 0.0f;
|
|
float dGamma, dBeta;
|
|
if(optimizer == 0)
|
|
{
|
|
dGamma = lr * gGamma + momentum * options[shift + BN_OPT_MG];
|
|
dBeta = lr * gBeta + momentum * options[shift + BN_OPT_MB];
|
|
options[shift + BN_OPT_MG] = dGamma;
|
|
options[shift + BN_OPT_MB] = dBeta;
|
|
}
|
|
else
|
|
{
|
|
const float mg = b1 * options[shift + BN_OPT_MG] + (1.0f - b1) * gGamma;
|
|
const float mb = b1 * options[shift + BN_OPT_MB] + (1.0f - b1) * gBeta;
|
|
const float vgOld = options[shift + BN_OPT_VG];
|
|
const float vbOld = options[shift + BN_OPT_VB];
|
|
const float vg = sqrt(b2 * vgOld * vgOld + (1.0f - b2) * gGamma * gGamma);
|
|
const float vb = sqrt(b2 * vbOld * vbOld + (1.0f - b2) * gBeta * gBeta);
|
|
dGamma = lt * mg / (vg > 0.0f ? vg : lt * 10.0f);
|
|
dBeta = lt * mb / (vb > 0.0f ? vb : lt * 10.0f);
|
|
options[shift + BN_OPT_MG] = mg;
|
|
options[shift + BN_OPT_MB] = mb;
|
|
options[shift + BN_OPT_VG] = vg;
|
|
options[shift + BN_OPT_VB] = vb;
|
|
}
|
|
dGamma = clamp(dGamma, -MAX_WEIGHT_DELTA, MAX_WEIGHT_DELTA);
|
|
dBeta = clamp(dBeta, -MAX_WEIGHT_DELTA, MAX_WEIGHT_DELTA);
|
|
// A non-finite delta means SKIP this unit's parameter write, never poison it - the moments above
|
|
// are already stored, matching the host's ordering exactly.
|
|
if(isnan(dGamma) || isinf(dGamma) || isnan(dBeta) || isinf(dBeta))
|
|
return;
|
|
options[shift + BN_OPT_GAMMA] = clamp(gamma + dGamma, -MAX_WEIGHT, MAX_WEIGHT);
|
|
options[shift + BN_OPT_BETA] = clamp(beta + dBeta, -MAX_WEIGHT, MAX_WEIGHT);
|
|
};
|
|
//+------------------------------------------------------------------+
|
|
//| LSTM (CNeuronLSTMOCL) - derived from scratch from the standard |
|
|
//| LSTM equations (NOT ported from any reference implementation). |
|
|
//| Single-timestep-truncated BPTT: gradient does not flow back into |
|
|
//| h_prev/c_prev from a previous step - only this step's gates and |
|
|
//| weights receive gradient. Adam-only optimizer. |
|
|
//| Weight layout: 4 gates (order: forget, input, output, candidate), |
|
|
//| each a block of H*(I+H+1) floats - per hidden unit: H weights on |
|
|
//| h_prev, I weights on the current input, 1 bias. |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void LSTM_Gates(__global float *matrix_w,
|
|
__global float *hidden_prev,
|
|
__global float *inputs,
|
|
__global float *concatenated,
|
|
int hidden_size, int input_size)
|
|
{
|
|
int id = get_global_id(0);
|
|
int gate = get_global_id(1);
|
|
int per_gate = hidden_size * (hidden_size + input_size + 1);
|
|
int shift = gate * per_gate + id * (hidden_size + input_size + 1);
|
|
float sum = 0.0f;
|
|
for(int k = 0; k < hidden_size; k++)
|
|
sum += hidden_prev[k] * matrix_w[shift + k];
|
|
for(int k = 0; k < input_size; k++)
|
|
sum += inputs[k] * matrix_w[shift + hidden_size + k];
|
|
sum += matrix_w[shift + hidden_size + input_size];
|
|
float val = (gate < 3) ? (1.0f / (1.0f + exp(-sum))) : tanh(sum);
|
|
concatenated[gate * hidden_size + id] = val;
|
|
}
|
|
//+------------------------------------------------------------------+
|
|
//| |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void LSTM_State(__global float *concatenated,
|
|
__global float *memory,
|
|
__global float *hidden_prev,
|
|
__global float *hidden_cache,
|
|
__global float *output,
|
|
int hidden_size)
|
|
{
|
|
int id = get_global_id(0);
|
|
float f = concatenated[id];
|
|
float ii = concatenated[hidden_size + id];
|
|
float o = concatenated[2 * hidden_size + id];
|
|
float g = concatenated[3 * hidden_size + id];
|
|
float c_prev = memory[id];
|
|
memory[hidden_size + id] = c_prev;
|
|
float c_t = f * c_prev + ii * g;
|
|
memory[id] = c_t;
|
|
hidden_cache[id] = hidden_prev[id];
|
|
output[id] = o * tanh(c_t);
|
|
}
|
|
//+------------------------------------------------------------------+
|
|
//| |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void LSTM_GateGradient(__global float *gradient,
|
|
__global float *memory,
|
|
__global float *concatenated,
|
|
__global float *concatenated_gradient,
|
|
int hidden_size)
|
|
{
|
|
int id = get_global_id(0);
|
|
float c_t = memory[id];
|
|
float c_prev = memory[hidden_size + id];
|
|
float f = concatenated[id];
|
|
float ii = concatenated[hidden_size + id];
|
|
float o = concatenated[2 * hidden_size + id];
|
|
float g = concatenated[3 * hidden_size + id];
|
|
float t = tanh(c_t);
|
|
float dh = gradient[id];
|
|
// Derivative floor on every gate below - see MIN_ACTIVATION_DERIVATIVE's declaration comment.
|
|
// LSTM gates are sigmoid/tanh by design and routinely sit near their saturated extremes (a fully
|
|
// "open" or "closed" gate), which is exactly where an unfloored derivative goes to zero and
|
|
// blocks gradient from flowing through the recurrent path at all.
|
|
float dc = dh * o * fmax(MIN_ACTIVATION_DERIVATIVE, 1 - t * t);
|
|
concatenated_gradient[2 * hidden_size + id] = dh * t * fmax(MIN_ACTIVATION_DERIVATIVE, o * (1 - o)); // output gate
|
|
concatenated_gradient[id] = dc * c_prev * fmax(MIN_ACTIVATION_DERIVATIVE, f * (1 - f)); // forget gate
|
|
concatenated_gradient[hidden_size + id] = dc * g * fmax(MIN_ACTIVATION_DERIVATIVE, ii * (1 - ii)); // input gate
|
|
concatenated_gradient[3 * hidden_size + id] = dc * ii * fmax(MIN_ACTIVATION_DERIVATIVE, 1 - g * g); // candidate
|
|
}
|
|
//+------------------------------------------------------------------+
|
|
//| |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void LSTM_WeightsGradient(__global float *concatenated_gradient,
|
|
__global float *hidden_cache,
|
|
__global float *inputs,
|
|
__global float *weights_gradient,
|
|
int hidden_size, int input_size)
|
|
{
|
|
int wi = get_global_id(0);
|
|
int per_gate = hidden_size * (hidden_size + input_size + 1);
|
|
int gate = wi / per_gate;
|
|
int rem = wi % per_gate;
|
|
int id = rem / (hidden_size + input_size + 1);
|
|
int k = rem % (hidden_size + input_size + 1);
|
|
float inp = (k < hidden_size ? hidden_cache[k] : (k < hidden_size + input_size ? inputs[k - hidden_size] : 1.0f));
|
|
weights_gradient[wi] = concatenated_gradient[gate * hidden_size + id] * inp;
|
|
}
|
|
//+------------------------------------------------------------------+
|
|
//| |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void LSTM_InputsGradient(__global float *concatenated_gradient,
|
|
__global float *matrix_w,
|
|
__global float *inputs_gradient,
|
|
int hidden_size, int input_size)
|
|
{
|
|
int j = get_global_id(0);
|
|
int per_gate = hidden_size * (hidden_size + input_size + 1);
|
|
float sum = 0.0f;
|
|
for(int gate = 0; gate < 4; gate++)
|
|
for(int id = 0; id < hidden_size; id++)
|
|
sum += concatenated_gradient[gate * hidden_size + id] * matrix_w[gate * per_gate + id * (hidden_size + input_size + 1) + hidden_size + j];
|
|
inputs_gradient[j] = sum;
|
|
}
|
|
//+------------------------------------------------------------------+
|
|
//| |
|
|
//+------------------------------------------------------------------+
|
|
//+------------------------------------------------------------------+
|
|
//| Sequence LSTM - ONE LAUNCH PER TIMESTEP. |
|
|
//| |
|
|
//| The recurrence is sequential: step t needs h/c from step t-1, and |
|
|
//| every hidden unit's gates read ALL of h_{t-1}. OpenCL can only |
|
|
//| barrier WITHIN a work-group, and nothing here controls how the |
|
|
//| runtime partitions the global size, so a single-launch loop with |
|
|
//| barrier() would be correct only by luck. Driving the loop from |
|
|
//| the host instead makes each launch an implicit global barrier - |
|
|
//| slower per step, correct on every device. |
|
|
//| |
|
|
//| Layouts match CPU_LSTMSeqForward in DirectML\WarriorCPU.cpp (the |
|
|
//| authoritative contract): gate order [f,i,o,g]; weight row `hid` |
|
|
//| of a gate is [H recurrent | Iw input | 1 bias]; caches are |
|
|
//| timestep-major. h_{-1} and c_{-1} are zero. |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void LSTM_SeqStepForward(__global float *matrix_w,
|
|
__global float *inputs,
|
|
__global float *cache_gates,
|
|
__global float *cache_cell,
|
|
__global float *cache_hidden,
|
|
__global float *output,
|
|
int hidden_size,
|
|
int step_inputs,
|
|
int steps,
|
|
int t)
|
|
{
|
|
int id = get_global_id(0);
|
|
int H = hidden_size;
|
|
int Iw = step_inputs;
|
|
int row = H + Iw + 1;
|
|
int per_gate = H * row;
|
|
__global float *x = inputs + t * Iw;
|
|
float g[4];
|
|
for(int gate = 0; gate < 4; gate++)
|
|
{
|
|
int shift = gate * per_gate + id * row;
|
|
float sum = 0;
|
|
if(t > 0)
|
|
{
|
|
__global float *h_prev = cache_hidden + (t - 1) * H;
|
|
for(int k = 0; k < H; k++)
|
|
sum += h_prev[k] * matrix_w[shift + k];
|
|
}
|
|
for(int k = 0; k < Iw; k++)
|
|
sum += x[k] * matrix_w[shift + H + k];
|
|
sum += matrix_w[shift + H + Iw];
|
|
g[gate] = (gate < 3) ? (1 / (1 + exp(-sum))) : tanh(sum);
|
|
cache_gates[t * 4 * H + gate * H + id] = g[gate];
|
|
}
|
|
float c_prev = (t > 0) ? cache_cell[(t - 1) * H + id] : 0;
|
|
float c_t = g[0] * c_prev + g[1] * g[3];
|
|
cache_cell[t * H + id] = c_t;
|
|
float h_t = g[2] * tanh(c_t);
|
|
cache_hidden[t * H + id] = h_t;
|
|
//--- the layer's visible output is the LAST hidden state
|
|
if(t == steps - 1)
|
|
output[id] = h_t;
|
|
}
|
|
//+------------------------------------------------------------------+
|
|
//| BPTT step 1/3: gate gradients for step t, and dc for step t-1. |
|
|
//| dh comes from the layer above at t == steps-1 and from the |
|
|
//| recurrence (dh_buf, written by LSTM_SeqStepInputGrad) otherwise. |
|
|
//| dc_buf must be zeroed by the host before the first step. |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void LSTM_SeqStepGateGrad(__global float *out_gradient,
|
|
__global float *dh_buf,
|
|
__global float *dc_buf,
|
|
__global float *cache_gates,
|
|
__global float *cache_cell,
|
|
__global float *gate_grad,
|
|
int hidden_size,
|
|
int steps,
|
|
int t)
|
|
{
|
|
int id = get_global_id(0);
|
|
int H = hidden_size;
|
|
__global float *gt = cache_gates + t * 4 * H;
|
|
float f = gt[id];
|
|
float ii = gt[H + id];
|
|
float o = gt[2 * H + id];
|
|
float g = gt[3 * H + id];
|
|
float c_t = cache_cell[t * H + id];
|
|
float c_prev = (t > 0) ? cache_cell[(t - 1) * H + id] : 0;
|
|
float tc = tanh(c_t);
|
|
float dh = (t == steps - 1) ? out_gradient[id] : dh_buf[id];
|
|
//--- the term the per-step kernels could never express: dc arriving from step t+1
|
|
float dcTot = dh * o * fmax(MIN_ACTIVATION_DERIVATIVE, 1 - tc * tc) + dc_buf[id];
|
|
gate_grad[2 * H + id] = dh * tc * fmax(MIN_ACTIVATION_DERIVATIVE, o * (1 - o));
|
|
gate_grad[id] = dcTot * c_prev * fmax(MIN_ACTIVATION_DERIVATIVE, f * (1 - f));
|
|
gate_grad[H + id] = dcTot * g * fmax(MIN_ACTIVATION_DERIVATIVE, ii * (1 - ii));
|
|
gate_grad[3 * H + id] = dcTot * ii * fmax(MIN_ACTIVATION_DERIVATIVE, 1 - g * g);
|
|
dc_buf[id] = dcTot * f;
|
|
}
|
|
//+------------------------------------------------------------------+
|
|
//| BPTT step 2/3: ACCUMULATE dW for step t. Accumulates rather than |
|
|
//| assigns because the gate weights are shared by every timestep - |
|
|
//| the host zeroes weights_gradient once before the backward loop. |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void LSTM_SeqStepWeightGrad(__global float *gate_grad,
|
|
__global float *cache_hidden,
|
|
__global float *inputs,
|
|
__global float *weights_gradient,
|
|
int hidden_size,
|
|
int step_inputs,
|
|
int t)
|
|
{
|
|
int wi = get_global_id(0);
|
|
int H = hidden_size;
|
|
int Iw = step_inputs;
|
|
int row = H + Iw + 1;
|
|
int per_gate = H * row;
|
|
int gate = wi / per_gate;
|
|
int rem = wi % per_gate;
|
|
int hid = rem / row;
|
|
int k = rem % row;
|
|
float inp;
|
|
if(k < H)
|
|
inp = (t > 0) ? cache_hidden[(t - 1) * H + k] : 0;
|
|
else
|
|
if(k < H + Iw)
|
|
inp = inputs[t * Iw + (k - H)];
|
|
else
|
|
inp = 1;
|
|
weights_gradient[wi] += gate_grad[gate * H + hid] * inp;
|
|
}
|
|
//+------------------------------------------------------------------+
|
|
//| BPTT step 3/3: dx_t (this step's slice of the input gradient) and |
|
|
//| dh_{t-1} (the recurrent path). One flat range over both so a |
|
|
//| narrow per-step input still saturates the device. |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void LSTM_SeqStepInputGrad(__global float *gate_grad,
|
|
__global float *matrix_w,
|
|
__global float *inputs_gradient,
|
|
__global float *dh_buf,
|
|
int hidden_size,
|
|
int step_inputs,
|
|
int t)
|
|
{
|
|
int j = get_global_id(0);
|
|
int H = hidden_size;
|
|
int Iw = step_inputs;
|
|
int row = H + Iw + 1;
|
|
int per_gate = H * row;
|
|
int col = (j < Iw) ? (H + j) : (j - Iw);
|
|
float sum = 0;
|
|
for(int gate = 0; gate < 4; gate++)
|
|
for(int hid = 0; hid < H; hid++)
|
|
sum += gate_grad[gate * H + hid] * matrix_w[gate * per_gate + hid * row + col];
|
|
if(j < Iw)
|
|
inputs_gradient[t * Iw + j] = sum;
|
|
else
|
|
dh_buf[j - Iw] = sum;
|
|
}
|
|
//+------------------------------------------------------------------+
|
|
//| |
|
|
//+------------------------------------------------------------------+
|
|
__kernel void LSTM_UpdateWeightsAdam(__global float *matrix_w,
|
|
__global float *weights_gradient,
|
|
__global float *matrix_m,
|
|
__global float *matrix_v,
|
|
float l, float b1, float b2)
|
|
{
|
|
int wi = get_global_id(0);
|
|
float g = weights_gradient[wi];
|
|
float mt = b1 * matrix_m[wi] + (1 - b1) * g;
|
|
// v squared back before the recursion - see UpdateWeightsAdam's comment above for why.
|
|
float vt = sqrt(b2 * matrix_v[wi] * matrix_v[wi] + (1 - b2) * g * g);
|
|
float delta = clamp(l * mt / (vt > 0 ? vt : l * 10) - l * WEIGHT_DECAY * matrix_w[wi], -MAX_WEIGHT_DELTA, MAX_WEIGHT_DELTA);
|
|
// Sign-agreement gate removed - see UpdateWeightsAdam's comment above for why.
|
|
matrix_w[wi] = clamp(matrix_w[wi] + delta, -MAX_WEIGHT, MAX_WEIGHT);
|
|
matrix_m[wi] = mt;
|
|
matrix_v[wi] = vt;
|
|
};
|
|
//+------------------------------------------------------------------+
|
|
//| SGD+momentum counterpart to LSTM_UpdateWeightsAdam above - same |
|
|
//| flat (already-elementwise) weights_gradient/matrix_dw layout, just|
|
|
//| the classic heavy-ball update instead of Adam's per-parameter step|
|
|
//+------------------------------------------------------------------+
|
|
__kernel void LSTM_UpdateWeightsMomentum(__global float *matrix_w,
|
|
__global float *weights_gradient,
|
|
__global float *matrix_dw,
|
|
float learning_rates, float momentum, int optimizer)
|
|
{
|
|
int wi = get_global_id(0);
|
|
// `optimizer` retained unused - see UpdateWeightsMomentum's comment above.
|
|
float delta = learning_rates * weights_gradient[wi] + momentum * matrix_dw[wi];
|
|
matrix_dw[wi] = delta;
|
|
matrix_w[wi] = clamp(matrix_w[wi] + delta, -MAX_WEIGHT, MAX_WEIGHT);
|
|
};
|