function D = create_demand(Gamma,GammaH,mass,massH,learning,learning_p)

if ~exist('learning','var')
    learning = @(t) 1;
    learning_p = @(t) 0;
end

D = struct(...
        'fn', @(r,rstar,q,qstar,t) 1/Gamma/learning(t) * (r - rstar) + mass*(r>rstar) - mass*(r<rstar),...
        'r', @(r,rstar,q,qstar,t) 1/Gamma/learning(t), ... % derivative w.r.t r
        'rstar', @(r,rstar,q,qstar,t) - 1/Gamma/learning(t), ... % derivative w.r.t rstar
        'q', @(r,rstar,q,qstar,t) 0, ... % derivative w.r.t. q
        'qstar', @(r,rstar,q,qstar,t) 0, ... % derivative w.r.t. qstar
        't', @(r,rstar,q,qstar,t) - 1/Gamma/learning(t)^2 * learning_p(t) * (r - rstar), ... % derivative w.r.t. t
        'explicit1',true,... % explicit solution for r to D + (r-r*)D_r = X ?
        'explicit1_fn',@(x,rstar,q,qstar,t) rstar + (x<-mass) * Gamma*learning(t)/2 * (x+mass) + (x > mass) * Gamma*learning(t)/2 * (x-mass),...
        ...'explicit2',true,... % explicit solution for r to D = X ?
        ...'explicit2_fn',@(x,rstar,q,qstar,t) rstar + (x<-mass) * Gamma * (x+mass) + (x > mass) * Gamma * (x-mass),...
        'plots', @(r,rstar,q,qstar,t) 1/Gamma/learning(t) * (r - rstar) + mass*(r>rstar) - mass*(r<rstar),... % function used to compute b_I along the transition
        'Hfn',@(r,rstar,q,qstar,t) 1/GammaH/learning(t) * (r - rstar) + massH*(r>rstar) - massH*(r<rstar),... % home intermediaries demand
        'Hr',@(r,rstar,q,qstar,t) 1/GammaH/learning(t) ,...
        'Hrstar',@(r,rstar,q,qstar,t) -1/GammaH/learning(t) ,...
        'Hq',@(r,rstar,q,qstar,t) 0 ,...
        'Hqstar',@(r,rstar,q,qstar,t) 0 ,...
        'Ht',@(r,rstar,q,qstar,t) - 1/GammaH/learning(t)^2 * learning_p(t) * (r - rstar) ,...
        'explicit3',true,... % explicit solution for r to D + DH = X ?
        'explicit3_fn',@(x,rstar,q,qstar,t) rstar + (x<-mass-massH) * learning(t)/(Gamma^(-1)+GammaH^(-1)) * (x+mass+massH) + (x > mass+massH) *  learning(t)/(Gamma^(-1)+GammaH^(-1)) * (x-mass-massH)...
        );
end