function [sawtooth_term,dsawtoothterm_diotagauss] = chi_sawtooth(q,q_squared,model,cp)

%% optional: sawteeth behaviour - only when chi_csawtooth is nonzero
% Construction of smooth F(q) function which increases transport in the
% core due to sawteeth activity

if cp.chi_csawtooth ~= 0
  sawteeth_exp  = exp( (q - 0.95)/cp.chi_wsawtooth); % q - 0.95 is done instead of q - 1, to influence the region q>1 less
  one_over_sawteeth_exp_plus_1_squared = (1+sawteeth_exp).^(-2);
  sawtooth_term = cp.chi_csawtooth ./ (1+sawteeth_exp); %F(q)
  dsawtoothterm_diotagauss = cp.chi_csawtooth * one_over_sawteeth_exp_plus_1_squared .* sawteeth_exp .* (cp.chi_wsawtooth.*q_squared);
else
  % must define on all execution paths
  sawtooth_term = ones(model.rgrid.nrhogauss,1);
  dsawtoothterm_diotagauss = zeros(model.rgrid.nrhogauss,1);
end

end