OpenEMS_MultiTermDebyeExperiments/AddDjordjevicSarkarMaterial.m at main · toammann/OpenEMS_MultiTermDebyeExperiments · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
function [CSX] = AddDjordjevicSarkarMaterial(varargin)
%
% Add a wideband dielectric material to a CSX struct using a multi-term
% Debye fit of the Djordjevic–Sarkar model.
%
% Calculates Debye parameters from a single (eps_r, tanD) measurement via
% 'calcDjordjevicSarkarApprox' and adds the material to the CSX struct.
%
% Input Parameters:
%   CSX          - CSX struct to which the material will be added
%   materialName - String name of the material
%
% Name-Value Parameters (passed to `calcDjordjevicSarkarApprox`):
%   'fMeas'      - Measurement frequency [Hz]
%   'epsRMeas'   - Relative permittivity ε_r at 'fMeas'
%   'tandMeas'   - Loss tangent tan(δ) at 'fMeas'
%   'f2'         - Upper corner frequency of the Djordjevic–Sarkar model [Hz]
%   'f1'         - Lower corner frequency [Hz] (if `lowFreqEvalType` = 0)
%   'epsRdc'     - Permittivity at DC (if `lowFreqEvalType` = 1)
%   'sigmaDC'    - Optional DC conductivity [S/m]
%   'nTermsPerDec' - Number of Debye terms per frequency decade
%
% Output:
%   CSX - Updated CSX struct including the defined wideband dielectric material
%
% Note:
%   - Internally uses 'calcDjordjevicSarkarApprox' to generate model parameters.
%   - See 'calcDjordjevicSarkarApprox' for detailed description of the model
%     fitting.
%
% Example:
%   CSX = AddDjordjevicSarkarMaterial(CSX, 'MyMaterial', ...
%     'fMeas', 1e9, 'epsRMeas', 4.2, 'tandMeas', 0.02, ...
%     'f1', 1e6, 'f2', 200e9);
%
% See also: calcDjordjevicSarkarApprox, AddDebyeMaterial
%
% Version History:
% v1.0  2025-06-30  Tobias Ammann  Initial version

  p = inputParser();
  p.FunctionName = 'AddDjordjevicSarkarMaterial';

  p.addRequired('CSX',      @isstruct);
  p.addRequired('matName',  @ischar);

  % Required Parameters
  % Set as Parameters, checked manually too match openEMS name-value pair style
  p.addParameter('fMeas',    [], @isPositiveScalar); % Hz
  p.addParameter('epsRMeas', [], @isPositiveScalar);
  p.addParameter('tandMeas', [], @isNonNegScalar);
  p.addParameter('f2',       [], @isPositiveScalar); % Upper fit frequency, Hz

  % Optional Parameters
  p.addParameter('lowFreqEvalType', 0,   @isIntegerScalar);   % 0 = use f1, 1 = use epsRdc
  p.addParameter('f1',              nan, @isPositiveScalar);  % lowFreqEvalType = 0, Hz
  p.addParameter('epsRdc',          inf, @isPositiveScalar);  % lowFreqEvalType = 1
  p.addParameter('sigmaDC',         0,   @isNonNegScalar);    % Siemens
  p.addParameter('nTermsPerDec',    1,   @(x) isScalar(x) && x >= 1); % Number of Debye terms per decade
  p.addParameter('plotEn',          0,   @isIntegerScalar);   % Enable/Disable plots of the model

  % Parse and manually verify required parameters
  p.parse(varargin{:});

  requiredParams = {'fMeas', 'epsRMeas', 'tandMeas', 'f2'};
  for i = 1:numel(requiredParams)
    param = requiredParams{i};
    if ismember(param, p.UsingDefaults)
      error('%s: Missing required parameter ''%s''.', p.FunctionName, param);
    end
  end

  if ((p.Results.lowFreqEvalType == 0) && ismember('f1', p.UsingDefaults))
    error(['%s: For ''lowFreqEvalType=0'' a value for f1 (Djordjevic Sarkar ',...
          'low corner frequency)must be specified.'], p.FunctionName);
  end

  if ((p.Results.lowFreqEvalType == 1) && ismember('epsRdc', p.UsingDefaults))
    error(['%s: For ''lowFreqEvalType=1'', a value for epsRdc (value of ',...
           'EpsilonR at DC) must be specified.'], p.FunctionName);
  end

   % Fit the model and receive Debye model parametes for openEMS
   paramDebye = calcDjordjevicSarkarApprox(...
      'fMeas', p.Results.fMeas,...
      'epsRMeas', p.Results.epsRMeas,...
      'tandMeas', p.Results.tandMeas,...
      'f1', p.Results.f1,...
      'f2', p.Results.f2,...
      'lowFreqEvalType', p.Results.lowFreqEvalType,...
      'epsRdc', p.Results.epsRdc,...
      'sigmaDC', p.Results.sigmaDC,...
      'nTermsPerDec', p.Results.nTermsPerDec,...
      'plotEn', p.Results.plotEn);

      CSX = p.Results.CSX;
      matName = p.Results.matName;

      CSX = AddDebyeMaterial(CSX, matName);
      CSX = SetMaterialProperty(CSX, matName, ...
                               'Epsilon', paramDebye.epsInf, ... % Epsilon here is acutally EpsilonRinf
                               'Kappa', paramDebye.sigmaDC);

      for i = 1:length(paramDebye.wi)

        CSX = SetMaterialProperty(CSX, matName,...
          ['EpsilonDelta_', int2str(i)],     paramDebye.deltaEpsT(i),...
          ['EpsilonRelaxTime_', int2str(i)], 1/paramDebye.wi(i));

      end
end

% Validation functions for input argument checks
function val = isNonNegScalar(x)
    val = isnumeric(x) && isscalar(x) && (x >= 0);
end

function val = isPositiveScalar(x)
    val = isnumeric(x) && isscalar(x) && (x > 0);
end

function val = isIntegerScalar(x)
    val = isnumeric(x) && isscalar(x) && (round(x) == x);
end

function val = isScalar(x)
    val = isnumeric(x) && isscalar(x);
end