diff --git a/docs/masks_and_checklist.txt b/docs/masks_and_checklist.txt new file mode 100644 index 0000000..7d3e3b8 --- /dev/null +++ b/docs/masks_and_checklist.txt @@ -0,0 +1,92 @@ +Simulink Block Masks and Completion Checklist + +Purpose: guide finishing the EMS_Solar_Rickshaw model. Each subsystem mask lists ports, key parameters, sensor mappings, and notes for hardware compatibility. + +1) PV_System (mask) +- Inports: Irradiance (W/m^2), Temperature (C) +- Outports: Vpv (V), Ipv (A), Ppv_kW, Vref (V) +- Internals: MATLAB Function calling src/pv_model and src/mppt_pando +- Notes: expose PV struct (module specs) as mask parameter; Vref is MPPT target voltage (to boost controller) + +2) Boost_Converter (mask) +- Inports: Vpv, Ipv, Vref +- Outports: Vdc_link, Ppv_out_kW +- Parameters: duty limits, switching freq (for power electronics block) +- Notes: implement converter loss model (efficiency), control loop to regulate Vpv->Vref (PI) or duty command from Vref + +3) DC_Link (mask) +- Inports: Ppv_in_kW, Pbidirectional_in_kW, Pload_dc_kW +- Outports: Vdc, P_net_kW +- Parameters: Cdc (capacitance), nominal Vdc +- Notes: model energy balance; provide Vdc sensor signal + +4) Battery (mask) +- Inports: P_batt_cmd_kW (positive discharge), Vbatt_measure +- Outports: SOC (0-1), Ibat (A), Vbatt (V) +- Parameters: nominal voltage, capacity_Ah, eff, SOC_init, SOC_min/max +- Internals: call src/battery_soc_update; expose coulomb-counting enable and initial SOC +- Hardware mapping: Vbatt from voltage divider, Ibat from bidirectional current sensor + +5) Bidirectional_Converter (mask) +- Inports: Vdc, Vbatt, Ibat_cmd (or P_cmd) +- Outports: P_batt_kW (positive discharge to DC bus), efficiency flags +- Parameters: charge/discharge efficiency, current limits +- Notes: control handles charge/discharge power; apply converter losses + +6) Grid (mask) +- Inports: Grid availability (bool), Grid command (import/export kW) +- Outports: P_grid_kW (positive = import), Vgrid (230V RMS), fault flag +- Notes: implement grid disconnect behavior; allow bidirectional flow + +7) Loads (mask) +- Inports: Residential_kW, Rickshaw_kW, rickshaw_allowed(bool) +- Outports: Total_load_kW, Critical_load_kW +- Notes: when rickshaw_allowed=false, set Rickshaw_kW=0; mark which loads are critical (residential critical) + +8) EMS (mask) +- Inports: Ppv_kW, Pload_kW, SOC, GridOk +- Outports: P_batt_cmd_kW, P_grid_import_kW, P_grid_export_kW, rickshaw_allowed, operating_mode (enum/string) +- Implementation: use src/ems_matlab_function; expose SOC limits and control hysteresis as mask parameters +- Hardware: operating_mode and commands map to MCU outputs (PWM duty, relay control) + +9) Measurement (mask) +- Inports: V/I sensors, Irradiance, Temp +- Outports: filtered signals, scaled sensor outputs, sample time parameter +- Notes: add anti-aliasing, simple low-pass filters, ADC scaling parameters + +10) Dashboard (mask) +- Inputs: key signals (Irradiance, Temp, Ppv, Vbatt, Ibat, SOC, Grid P, Loads, Mode) +- Implement: Scopes, To Workspace blocks, and display annotations + +Completion Checklist +1. Run scripts/init_simulation in MATLAB to populate PV and BATT structs. +2. Run scripts/import_nasa_power('full\path\to\CSV') to create timetable and timeseries. +3. Run scripts/generate_load_profiles to get 24h loads. +4. Run scripts/build_simulink_model then scripts/flesh_out_model to create skeleton and populate function blocks. +5. Open EMS_Solar_Rickshaw in Simulink. Inspect each subsystem created by the script. +6. For each subsystem: + - Add PI controllers (Boost & Bidirectional) and tune gains for 1-hour step resolution. + - Implement converter efficiency/loss models (set 95% typical values). + - Add measurement filters and ADC scaling blocks. + - Wire EMS outputs to converter control inputs and grid relays. +7. Replace placeholder MATLAB Function blocks where necessary with S-Functions or Simscape Power Systems blocks if desired for more fidelity. +8. Validate energy balance: PV generation = load + battery charge + grid export (consider losses). +9. Test scenarios: + - Clear-sky day (use provided NASA file) — expect exports when SOC hits max. + - Cloudy day (low irradiance) — battery discharges, grid imports when SOC hits min. + - Grid failure (set GridOk=false) — rickshaw_allowed=false, battery supplies critical load only. + - Hysteresis test: vary SOC near limits to ensure stable switching (add small deadband/hysteresis parameters). +10. Generate plots: irradiance vs time, PV power, SOC, battery power, grid import/export, loads, operating mode timeline. + +Tuning & Validation Tips +- Use sample time = 3600s for high-level energy simulation. For converter control dynamics, consider a separate faster-timescale model if needed. +- Keep units consistent (kW, V, A, Ah). Document signal units in block masks. +- Expose key parameters as mask variables for repeatable experiments. + +Deliverables to produce after finishing model +- Simulink model file (.slx) +- MATLAB init and import scripts (already created) +- EMS MATLAB Function block source (src/ems_matlab_function.m) +- Plots and short report describing scenarios and results + +If desired, next action: create mask parameter templates (Simulink mask XML snippets) or commit these files and open a PR. Which is preferred? \ No newline at end of file diff --git a/docs/signal_mapping.txt b/docs/signal_mapping.txt new file mode 100644 index 0000000..88f5734 --- /dev/null +++ b/docs/signal_mapping.txt @@ -0,0 +1,17 @@ +Signal mapping (for future hardware sensors) + +- PV Irradiance: sensor -> pyranometer (raw W/m^2) -> analog input +- PV Temperature: NTC/thermistor on module -> analog input +- PV Voltage: voltage divider -> ADC -> Vpv (Volts) +- PV Current: shunt or Hall sensor -> Ipv (Amps) +- Battery Voltage: voltage divider -> ADC -> Vbatt (Volts) +- Battery Current: bidirectional current sensor (Hall/shunt) -> Ibat (Amps) +- Battery SOC: computed in MCU from Coulomb counting + voltage lookup (signal: SOC, 0-1) +- Grid Voltage: voltage transformer -> ADC -> Vgrid (Volts) +- Grid Current: current transformer -> ADC -> Igrid (Amps) + +Which signals originate from sensors vs estimated: +- Sensors: Vpv, Ipv, Vbatt, Ibat, Vgrid, Igrid, Irradiance, ModuleTemp +- Estimated: SOC (from MCU algorithm) + +Recommended naming in Simulink: prefix signals with src_ (sensor), est_ (estimated), ctrl_ (control commands) diff --git a/scripts/add_masks_and_controllers.m b/scripts/add_masks_and_controllers.m new file mode 100644 index 0000000..f4fe1fe --- /dev/null +++ b/scripts/add_masks_and_controllers.m @@ -0,0 +1,113 @@ +% add_masks_and_controllers Create masks for subsystems and add simple PI controllers +% Run in MATLAB with Simulink available after building and fleshing out the model. + +modelName = 'EMS_Solar_Rickshaw'; +if ~bdIsLoaded(modelName) + try + open_system(modelName); + catch + error('Model %s not found. Run scripts/build_simulink_model and scripts/flesh_out_model first.', modelName); + end +end + +% Helper to safely create mask and parameters +function create_simple_mask(sysPath, params) + try + maskObj = Simulink.Mask.create(sysPath); + catch + maskObj = Simulink.Mask.get(sysPath); + end + for k=1:numel(params) + p = params{k}; + try + Simulink.Mask.addParameter(maskObj,p{:}); + catch + % ignore duplicates + end + end +end + +% Define masks and parameters +create_simple_mask([modelName '/PV_System'], { + {'edit','Module_P_rated','445'}, + {'edit','Module_Vmp','33.64'}, + {'edit','Module_Voc','39.65'} +}); + +create_simple_mask([modelName '/Boost_Converter'], { + {'edit','Boost_Kp','0.1'}, + {'edit','Boost_Ki','0.01'}, + {'edit','Boost_Duty_Max','0.95'}, + {'edit','Boost_Duty_Min','0.05'} +}); + +create_simple_mask([modelName '/Bidirectional_Converter'], { + {'edit','Bi_Kp','0.2'}, + {'edit','Bi_Ki','0.02'}, + {'edit','Bi_Imax_A','100'} +}); + +create_simple_mask([modelName '/Battery'], { + {'edit','Batt_Vnom','72'}, + {'edit','Batt_Capacity_Ah','100'}, + {'edit','Batt_SOC_init','0.6'} +}); + +create_simple_mask([modelName '/Grid'], { + {'edit','Grid_Vrms','230'}, + {'edit','Grid_freq','50'} +}); + +create_simple_mask([modelName '/EMS'], { + {'edit','SOC_min','0.20'}, + {'edit','SOC_max','0.95'}, + {'edit','Hysteresis_pct','0.02'} +}); + +% Add a simple PI controller to Boost_Converter: Duty = Kp*(Vref-Vpv) + Ki*integral(Vref-Vpv) +boostPath = [modelName '/Boost_Converter']; +try + delete_block([boostPath '/*']); +catch +end +add_block('built-in/Inport',[boostPath '/In_Vpv'],'Position',[30 30 60 50]); +add_block('built-in/Inport',[boostPath '/In_Vref'],'Position',[30 80 60 100]); +add_block('built-in/Sum',[boostPath '/SumErr'],'Position',[110 40 140 80],'Inputs','+-'); +add_block('built-in/Gain',[boostPath '/Kp'],'Position',[200 30 230 60],'Gain','Boost_Kp'); +add_block('built-in/Integrator',[boostPath '/Int'],'Position',[200 80 230 120]); +add_block('built-in/Gain',[boostPath '/Ki'],'Position',[300 80 330 120],'Gain','Boost_Ki'); +add_block('built-in/Sum',[boostPath '/SumOut'],'Position',[360 50 390 100],'Inputs','++'); +add_block('built-in/Saturation',[boostPath '/Sat'],'Position',[440 50 470 90],'UpperLimit','Boost_Duty_Max','LowerLimit','Boost_Duty_Min'); +add_block('built-in/Outport',[boostPath '/Out_Duty'],'Position',[520 60 550 80]); + +add_line(boostPath,'In_Vref/1','SumErr/1'); +add_line(boostPath,'In_Vpv/1','SumErr/2'); +add_line(boostPath,'SumErr/1','Kp/1'); +add_line(boostPath,'SumErr/1','Int/1'); +add_line(boostPath,'Int/1','Ki/1'); +add_line(boostPath,'Kp/1','SumOut/1'); +add_line(boostPath,'Ki/1','SumOut/2'); +add_line(boostPath,'SumOut/1','Sat/1'); +add_line(boostPath,'Sat/1','Out_Duty/1'); + +% Add simple controller to Bidirectional_Converter to follow requested battery power (P_cmd) +biPath = [modelName '/Bidirectional_Converter']; +try + delete_block([biPath '/*']); +catch +end +add_block('built-in/Inport',[biPath '/In_Vdc'],'Position',[30 30 60 50]); +add_block('built-in/Inport',[biPath '/In_Vbatt'],'Position',[30 80 60 100]); +add_block('built-in/Inport',[biPath '/In_Pcmd_kW'],'Position',[30 130 60 150]); +add_block('built-in/Gain',[biPath '/PtoI'],'Position',[200 120 240 160],'Gain','1/(Batt_Vnom)'); +add_block('built-in/Saturation',[biPath '/Ilimit'],'Position',[300 120 340 160],'UpperLimit','Bi_Imax_A','LowerLimit','-Bi_Imax_A'); +add_block('built-in/Outport',[biPath '/Out_Ibat'], 'Position',[420 130 450 150]); +add_line(biPath,'In_Pcmd_kW/1','PtoI/1'); +add_line(biPath,'PtoI/1','Ilimit/1'); +add_line(biPath,'Ilimit/1','Out_Ibat/1'); + +% Save model +save_system(modelName); +fprintf('Masks and basic controllers added. Tune mask parameters in Simulink UI.\n'); + +end diff --git a/scripts/build_simulink_model.m b/scripts/build_simulink_model.m new file mode 100644 index 0000000..84500fb --- /dev/null +++ b/scripts/build_simulink_model.m @@ -0,0 +1,30 @@ +% build_simulink_model Programmatically create a Simulink model skeleton for the EMS project +% Run this from MATLAB with Simulink available. It creates subsystems for the required components. + +modelName = 'EMS_Solar_Rickshaw'; +if bdIsLoaded(modelName) + close_system(modelName,0); +end + +new_system(modelName); +open_system(modelName); + +% Helper to add a subsystem +addSub = @(name, pos) add_block('built-in/Subsystem', [modelName '/' name], 'Position', pos); + +% Add subsystems with approximate positions +addSub('PV_System',[30 30 180 140]); +addSub('Boost_Converter',[220 30 370 140]); +addSub('DC_Link',[410 30 540 140]); +addSub('Battery',[30 200 180 320]); +addSub('Bidirectional_Converter',[220 200 370 320]); +addSub('Grid',[600 30 740 140]); +addSub('Loads',[600 200 740 320]); +addSub('EMS',[410 200 540 320]); +addSub('Measurement',[410 350 540 440]); +addSub('Dashboard',[600 350 900 600]); + +% Save model +save_system(modelName); +fprintf('Created Simulink model skeleton: %s\n', modelName); +end \ No newline at end of file diff --git a/scripts/flesh_out_model.m b/scripts/flesh_out_model.m new file mode 100644 index 0000000..78f103e --- /dev/null +++ b/scripts/flesh_out_model.m @@ -0,0 +1,74 @@ +% flesh_out_model Populate the EMS_Solar_Rickshaw Simulink model with functional blocks +% Run this in MATLAB (Simulink must be available). The script adds MATLAB Function blocks +% that call src/pv_model, src/mppt_pando, src/battery_soc_update and hooks them into the skeleton. + +modelName = 'EMS_Solar_Rickshaw'; +if ~bdIsLoaded(modelName) + try + open_system(modelName); + catch + error('Model %s not found. Run scripts/build_simulink_model.m first.', modelName); + end +end + +% Add PV_System internals: In ports for Irradiance, Temperature, Out ports for Vpv, Ipv, Ppv +pvPath = [modelName '/PV_System']; +% Remove existing contents for idempotence +delete_block([pvPath '/*']); +add_block('built-in/Inport', [pvPath '/In_Irradiance'], 'Position',[30 30 60 50], 'Port', '1'); +add_block('built-in/Inport', [pvPath '/In_Temp'], 'Position',[30 80 60 100], 'Port', '2'); +mf = add_block('built-in/MATLABFcn', [pvPath '/pv_model_func'], 'Position',[130 30 270 120]); +set_param(mf, 'FunctionName', 'pv_model'); +add_block('built-in/Outport', [pvPath '/Out_Vpv'], 'Position',[320 20 350 40], 'Port', '1'); +add_block('built-in/Outport', [pvPath '/Out_Ipv'], 'Position',[320 70 350 90], 'Port', '2'); +add_block('built-in/Outport', [pvPath '/Out_Ppv_kW'], 'Position',[320 120 350 140], 'Port', '3'); + +% Connect lines +add_line(pvPath, 'In_Irradiance/1','pv_model_func/1'); +add_line(pvPath, 'In_Temp/1','pv_model_func/2'); +% Note: pv_model requires PV struct in base workspace. The MATLAB Function block will call: [Vmp,Imp,Ppv,Vpv,Ipv] = pv_model(u(1),u(2),PV); + +% Add MPPT block inside PV_System +mpptBlk = add_block('built-in/MATLABFcn', [pvPath '/mppt_pando'], 'Position',[130 140 270 220]); +set_param(mpptBlk, 'FunctionName', 'mppt_pando'); +add_block('built-in/Outport', [pvPath '/Out_Vref'], 'Position',[320 160 350 180], 'Port', '4'); + +% Loads subsystem: attach Inports for residential and rickshaw loads and Outport for total load +loadPath = [modelName '/Loads']; +delete_block([loadPath '/*']); +add_block('built-in/Inport', [loadPath '/In_res_kW'], 'Position',[30 30 60 50]); +add_block('built-in/Inport', [loadPath '/In_rick_kW'], 'Position',[30 80 60 100]); +add_block('built-in/MathFunction', [loadPath '/SumLoads'], 'Position',[130 40 200 120], 'Operator', 'sum'); +add_block('built-in/Outport', [loadPath '/Out_total_kW'], 'Position',[320 60 350 80]); +add_line(loadPath, 'In_res_kW/1','SumLoads/1'); +add_line(loadPath, 'In_rick_kW/1','SumLoads/2'); +add_line(loadPath, 'SumLoads/1','Out_total_kW/1'); + +% Battery subsystem: simple interface for SOC update +batPath = [modelName '/Battery']; +delete_block([batPath '/*']); +add_block('built-in/Inport', [batPath '/In_Pcmd_kW'], 'Position',[30 30 60 50]); +add_block('built-in/Inport', [batPath '/In_Vbatt'], 'Position',[30 80 60 100]); +batFn = add_block('built-in/MATLABFcn', [batPath '/battery_soc_update'], 'Position',[130 30 300 140]); +set_param(batFn, 'FunctionName', 'battery_soc_update'); +add_block('built-in/Outport', [batPath '/Out_SOC'], 'Position',[350 40 380 60]); +add_block('built-in/Outport', [batPath '/Out_Ibat'], 'Position',[350 90 380 110]); + +% EMS subsystem: replace with MATLAB Function call to src/ems_matlab_function.m +emsPath = [modelName '/EMS']; +delete_block([emsPath '/*']); +add_block('built-in/Inport', [emsPath '/In_PV_kW'], 'Position',[30 30 60 50]); +add_block('built-in/Inport', [emsPath '/In_Load_kW'], 'Position',[30 80 60 100]); +add_block('built-in/Inport', [emsPath '/In_SOC'], 'Position',[30 130 60 150]); +add_block('built-in/Inport', [emsPath '/In_GridOk'], 'Position',[30 180 60 200]); +emsFn = add_block('built-in/MATLABFcn', [emsPath '/ems_function'], 'Position',[130 30 380 200]); +set_param(emsFn, 'FunctionName', 'ems_matlab_function'); +add_block('built-in/Outport', [emsPath '/Out_grid_import_kW'], 'Position',[420 30 450 50]); +add_block('built-in/Outport', [emsPath '/Out_grid_export_kW'], 'Position',[420 80 450 100]); +add_block('built-in/Outport', [emsPath '/Out_batt_power_kW'], 'Position',[420 130 450 150]); +add_block('built-in/Outport', [emsPath '/Out_rickshaw_allowed'], 'Position',[420 180 450 200]); + +% Save model +save_system(modelName); +fprintf('Fleshed out model subsystems in %s. Open in Simulink to complete wiring and tune parameters.\n', modelName); +end \ No newline at end of file diff --git a/scripts/generate_load_profiles.m b/scripts/generate_load_profiles.m new file mode 100644 index 0000000..76a0db7 --- /dev/null +++ b/scripts/generate_load_profiles.m @@ -0,0 +1,41 @@ +function loads = generate_load_profiles(startDatetime) +% generate_load_profiles Create 24-hour realistic residential and rickshaw load profiles +% loads = generate_load_profiles(startDatetime) +% startDatetime (datetime) optional, default uses today at 00:00 + +if nargin<1 || isempty(startDatetime) + startDatetime = dateshift(datetime('now'),'start','day'); +end + +hours = startDatetime + hours(0:23); + +% Residential load (kW) - base + morning and evening peaks +base = 0.6; % kW base +morning_peak = [0.6 0.8 1.0 1.3 1.5 1.4 1.2 1.0]; % early morning ramp +evening_peak = [1.2 1.4 1.6 1.8 1.6 1.3 1.0 0.8]; % evening +res = zeros(1,24) + base; +res(6:13) = res(6:13) + [0.2 0.4 0.6 0.6 0.5 0.4 0.3 0.2]; % day usage +res(18:23) = res(18:23) + [0.6 0.8 1.2 1.4 1.2 1.0]; +res = res; % kW + +% Electric rickshaw charging profile (kW) - assume multiple rickshaws, charging mostly daytime and afternoon +rickshaw = zeros(1,24); +rickshaw(8:11) = 0.5; % morning quick charging demand +rickshaw(13:16) = 0.8; % midday charging +rickshaw(17:19) = 0.6; % evening top-ups + +% Add some randomness to simulate variability but keep reproducible +rng(0); +res = res .* (1 + 0.05*randn(size(res))); +rickshaw = rickshaw .* (1 + 0.05*randn(size(rickshaw))); + +% Ensure non-negative +res(res<0) = 0; +rickshaw(rickshaw<0) = 0; + +% Create timetable +loads.t = hours'; +loads.residential_kW = res'; +loads.rickshaw_kW = rickshaw'; +loads.timetable = timetable(loads.t, loads.residential_kW, loads.rickshaw_kW, 'VariableNames', {'Residential_kW','Rickshaw_kW'}); +end \ No newline at end of file diff --git a/scripts/import_nasa_power.m b/scripts/import_nasa_power.m new file mode 100644 index 0000000..a5e27dc --- /dev/null +++ b/scripts/import_nasa_power.m @@ -0,0 +1,86 @@ +function data = import_nasa_power(csvPath) +% import_nasa_power Read NASA POWER CSV and return timetable and timeseries +% data = import_nasa_power(csvPath) +% Expects CSV with header block and a header line starting with YEAR,MO,DY,HR,... + +if nargin<1 || isempty(csvPath) + [file, path] = uigetfile({'*.csv'}, 'Select NASA POWER CSV'); + if isequal(file,0), error('No file selected'); end + csvPath = fullfile(path,file); +end + +fid = fopen(csvPath,'r'); +if fid==-1, error('Could not open %s', csvPath); end + +% Find header line containing column names +line = ''; +lineNum = 0; +while ischar(line) + line = fgetl(fid); + lineNum = lineNum + 1; + if ~ischar(line), break; end + if startsWith(strtrim(line),'YEAR,') + headerLine = line; + break; + end +end + +if ~exist('headerLine','var') + fclose(fid); + error('Could not find data header line in CSV'); +end + +% Rewind and read from the header line using readtable with NumHeaderLines +fclose(fid); +numHeaderLines = lineNum - 1; +opts = detectImportOptions(csvPath,'NumHeaderLines',numHeaderLines);% may infer types +T = readtable(csvPath, opts); + +% Build datetime vector +years = T.YEAR; months = T.MO; days = T.DY; hours = T.HR; +% NASA LST hours are 0-23; create datetime in local timezone (no tz attached) +dt = datetime(years,months,days,hours,0,0); + +% Extract variables (column names might vary slightly) +varNames = T.Properties.VariableNames; +% Look for irradiance and temperature columns +irrName = intersect(varNames, {'ALLSKY_SFC_SW_DWN','ALLSKY_SFC_SW_DWN..W.m.2.'}); +tempName = intersect(varNames, {'T2M','T2M..C.'}); +if isempty(irrName) || isempty(tempName) + % try by position (last two columns) + irr = T{:,end-1}; + temp = T{:,end}; +else + irr = T{:,irrName}; + temp = T{:,tempName}; +end + +% Convert irradiance units: NASA gives Wh/m^2 per hour; convert to W/m^2 by dividing by 1 hour (3600s) +% But hourly total irradiance in Wh/m^2 -> average W/m^2 over the hour = Wh/m^2 +% For simulation at 1-hour steps, use Wh/m^2 as proxy or convert to W/m^2 by leaving as-is as representative irradiance +% Here keep as W/m^2 (approx) by leaving numeric value (common practice uses W/m^2 approx equal to Wh/m^2 for hourly average) +irradiance = double(irr); +temperature = double(temp); + +% Create timetable +TT = timetable(dt, irradiance, temperature,'VariableNames',{'Irradiance','Temperature'}); + +% Create timeseries objects for convenience +tsIrr = timeseries(TT.Irradiance, TT.dt, 'Name', 'Irradiance'); +tsIrr.TimeInfo.Units = 'hours'; +tsTemp = timeseries(TT.Temperature, TT.dt, 'Name', 'Temperature'); +tsTemp.TimeInfo.Units = 'hours'; + +data.table = T; +data.timetable = TT; +data.tsIrradiance = tsIrr; +data.tsTemperature = tsTemp; + +% Save a MAT for quick reload +[folder, name] = fileparts(csvPath); +try + save(fullfile(folder,[name '.mat']), 'TT', 'T'); +catch + % ignore save errors +end +end \ No newline at end of file diff --git a/scripts/init_simulation.m b/scripts/init_simulation.m new file mode 100644 index 0000000..025ba1a --- /dev/null +++ b/scripts/init_simulation.m @@ -0,0 +1,44 @@ +% init_simulation Initialize parameters and load data for EMS Solar-Rickshaw simulation +% Run this before opening/running the Simulink model + +% Simulation timing +Ts = 3600; % sample time: 1 hour in seconds +simHours = 24; +stopTime = num2str(simHours*3600); + +% PV module specs (JA Solar JAM54D41 MB 445W) +PV.module.P_rated = 445; % W +PV.module.Voc = 39.65; % V +PV.module.Isc = 13.90; % A +PV.module.Vmp = 33.64; % V +PV.module.Imp = 13.23; % A +PV.module.eta = 0.228; % module efficiency +PV.module.alpha_P = -0.00290; % Pmax coeff per degC (fraction per degC) +PV.module.beta_Voc = -0.00250; % Voc coeff per degC +PV.module.alpha_Isc = 0.00045; % Isc coeff per degC + +% Battery specs +BATT.nominal_voltage = 72; % V +BATT.capacity_Ah = 100; % choose capacity so energy ~7.2 kWh (example) -> settable by user +BATT.initial_SOC = 0.60; % 60% +BATT.SOC_min = 0.20; +BATT.SOC_max = 0.95; +BATT.efficiency = 0.95; % round-trip + +% Grid specs +GRID.Vrms = 230; GRID.freq = 50; GRID.available = true; + +% Assign to base workspace for Simulink +assignin('base','Ts',Ts); +assignin('base','stopTime',stopTime); +assignin('base','PV',PV); +assignin('base','BATT',BATT); +assignin('base','GRID',GRID); + +% Load NASA POWER data (prompt user to select CSV if not present) +% Example usage: data = import_nasa_power('path_to_csv.csv'); + +disp('Initialization complete. Run import_nasa_power and generate_load_profiles as needed.'); + +clear Ts stopTime PV BATT GRID; % keep in base via assignin +end \ No newline at end of file diff --git a/src/battery_soc_update.m b/src/battery_soc_update.m new file mode 100644 index 0000000..f66aff1 --- /dev/null +++ b/src/battery_soc_update.m @@ -0,0 +1,32 @@ +function [SOC_new, Ibat_A] = battery_soc_update(SOC_prev, P_batt_kW, Vbatt, dt_sec, BATT) +% battery_soc_update Update battery SOC using simple Coulomb-counting +% Inputs: +% SOC_prev (0-1), P_batt_kW (positive = discharge to load, negative = charge), +% Vbatt (V), dt_sec (seconds), BATT struct with fields capacity_Ah and efficiency +% Outputs: SOC_new (0-1), Ibat_A (signed, positive = discharge) + +% Protect inputs +if nargin<5 || isempty(BATT) + error('BATT struct with capacity_Ah and efficiency required'); +end + +% Compute current (A). Add small eps to avoid div by zero +I = (P_batt_kW*1000)/(Vbatt + eps); +% I positive -> discharge, I negative -> charging +Ah_delta = I * (dt_sec/3600); + +% Account for efficiency: when charging, store less energy; when discharging, available energy reduced +if I < 0 + % charging: effective stored Ah reduced by efficiency + Ah_stored = -Ah_delta * BATT.efficiency; + SOC_new = SOC_prev + (Ah_stored / BATT.capacity_Ah); +else + % discharging + Ah_used = Ah_delta / BATT.efficiency; + SOC_new = SOC_prev - (Ah_used / BATT.capacity_Ah); +end + +% Clamp +SOC_new = min(max(SOC_new,0),1); +Ibat_A = I; +end \ No newline at end of file diff --git a/src/ems_matlab_function.m b/src/ems_matlab_function.m new file mode 100644 index 0000000..4bb92f1 --- /dev/null +++ b/src/ems_matlab_function.m @@ -0,0 +1,60 @@ +function [mode, grid_import_kW, grid_export_kW, batt_power_kW, rickshaw_allowed] = ems_matlab_function(pv_power_kW, load_power_kW, soc, grid_ok) +% EMS MATLAB Function implementing rule-based energy management +% Inputs: pv_power_kW, load_power_kW, soc (0-1), grid_ok (boolean) +% Outputs: mode (string), grid_import_kW, grid_export_kW, batt_power_kW (positive = discharge to load), rickshaw_allowed (boolean) + +% Initialize outputs +mode = "Undefined"; +grid_import_kW = 0; +grid_export_kW = 0; +batt_power_kW = 0; % positive means battery -> load +rickshaw_allowed = true; + +% Safety bounds +SOC_MIN = 0.20; +SOC_MAX = 0.95; + +if ~grid_ok + % Case 5: Grid failure -> island mode + mode = "Island"; + rickshaw_allowed = false; % disable rickshaw charging + % Battery supplies critical load only (assume load_power_kW is critical) + if soc > SOC_MIN + % discharge battery to supply load + batt_power_kW = min(load_power_kW, (soc - SOC_MIN)*1000); % placeholder scaling; final model uses Ah and V + else + % battery depleted, cannot supply + grid_import_kW = 0; batt_power_kW = 0; + end + return; +end + +% Normal operation +if pv_power_kW > load_power_kW + excess = pv_power_kW - load_power_kW; + if soc < SOC_MAX + % Case 1: charge battery with excess + mode = "PV>Load: Charge Battery"; + % Send excess to battery (after converter losses handled elsewhere) + batt_power_kW = -excess; % negative means charging + grid_export_kW = 0; + else + % Case 2: export to grid + mode = "PV>Load: Export to Grid"; + grid_export_kW = excess; + batt_power_kW = 0; + end +else + deficit = load_power_kW - pv_power_kW; + if soc > SOC_MIN + % Case 3: battery supplies load + mode = "PV prevP_W + % if power increased, keep direction + % no change +else + % power decreased -> reverse direction + dir = -dir; +end + +Vref = Vmeas + dir*step; +if Vref < 0, Vref = 0; end + +% update stored values +prevP = Pmeas; +prevV = Vmeas; + +status.step = step; +status.direction = dir; +end \ No newline at end of file diff --git a/src/pv_model.m b/src/pv_model.m new file mode 100644 index 0000000..47c58db --- /dev/null +++ b/src/pv_model.m @@ -0,0 +1,31 @@ +function [Vmp_est, Imp_est, Ppv_kW, Vpv, Ipv] = pv_model(irradiance_Wpm2, tempC, PV) +% pv_model Simple performance model for a PV module (educational) +% Inputs: irradiance_Wpm2, tempC, PV struct (from init_simulation) +% Outputs: Vmp_est (V), Imp_est (A), Ppv_kW, Vpv (V), Ipv (A) + +% Guard inputs +G = max(irradiance_Wpm2,0); +T = tempC; + +% Use module efficiency and P_rated to estimate area +eta = PV.module.eta; +P_rated = PV.module.P_rated; % W at STC +area_m2 = P_rated/(eta*1000 + eps); + +% Estimated module power at given irradiance and temp (simple linear temp coeff) +P_module_W = P_rated * (G/1000) * (1 + PV.module.alpha_P*(T-25)); +P_module_W = max(P_module_W,0); + +% Estimate Vmp change with temperature using Voc coeff approx +Vmp_est = PV.module.Vmp * (1 + PV.module.beta_Voc*(T-25)); +if Vmp_est <= 0, Vmp_est = PV.module.Vmp; end + +% Estimate Imp from P and Vmp +Imp_est = P_module_W / (Vmp_est + eps); + +% Provide instantaneous PV voltage/current approx +Vpv = Vmp_est; +Ipv = Imp_est; + +Ppv_kW = P_module_W/1000; +end \ No newline at end of file